MB1 CDI Follow-up Confirmatory Analysis
The ManyBabies Analysis Team
Wed Jul 21 16:01:53 2021
1 Introduction
In this script, we go through all the pre-registered proposed analyses. As a reminder, the main research questions where as follows:
- To what extent does infants#’ preference for IDS as measured in a laboratory setting predict their vocabulary at 18 and 24 months?
- Does the relation between IDS preference and vocabulary size change over development?
- Are there systematic differences in the strength of this relationship across the language communities in our sample?
We first present the main “sample theory based” analyses (also known as frequentist), separately on the North American and UK samples in parallel to answer our first two research questions, then on the total dataset to answer our third research question. We then provide additional Bayesian statistics where a null effect was found, as specified in the pre-registration.
# Library imports, general settings ==============
library(tidyverse); library(egg)
library(knitr)
library(lme4); library(lmerTest); library(simr)
# As in our discussion with Mike, we will use lmerTest for calculating p values
# in the mixed-effects models (noted as a deviation)
library(brms); library(rstan)
library(future)
plan(multicore, workers = 8) # Swap to multiprocess if running from Rstudio
library(lattice)
library(effects)
library(sjPlot)
library(robustlmm)
library(car)
# Load model comparison functions
source("helper/lrtests.R")
# Deal with package priority issues
select <- dplyr::select
theme_set(theme_bw(base_size = 10))
options("future" = T)
#knitr::opts_chunk$set(cache = TRUE)
print(sessionInfo()) #listing all info about R and packages info## R version 4.0.3 (2020-10-10)
## Platform: x86_64-apple-darwin17.0 (64-bit)
## Running under: macOS Catalina 10.15.7
##
## Matrix products: default
## BLAS: /Library/Frameworks/R.framework/Versions/4.0/Resources/lib/libRblas.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.0/Resources/lib/libRlapack.dylib
##
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] furrr_0.2.1 future.apply_1.7.0 bridgesampling_1.0-0
## [4] car_3.0-10 robustlmm_2.4-4 sjPlot_2.8.7
## [7] effects_4.2-0 carData_3.0-4 lattice_0.20-41
## [10] future_1.21.0 rstan_2.21.2 StanHeaders_2.21.0-7
## [13] brms_2.14.4 Rcpp_1.0.5 simr_1.0.5
## [16] lmerTest_3.1-3 lme4_1.1-26 Matrix_1.3-2
## [19] knitr_1.30 egg_0.4.5 gridExtra_2.3
## [22] forcats_0.5.0 stringr_1.4.0 dplyr_1.0.2
## [25] purrr_0.3.4 readr_1.4.0 tidyr_1.1.2
## [28] tibble_3.0.4 ggplot2_3.3.3 tidyverse_1.3.0
##
## loaded via a namespace (and not attached):
## [1] tidyselect_1.1.0 htmlwidgets_1.5.3 grid_4.0.3
## [4] munsell_0.5.0 codetools_0.2-18 effectsize_0.4.1
## [7] statmod_1.4.35 DT_0.17 miniUI_0.1.1.1
## [10] withr_2.3.0 Brobdingnag_1.2-6 colorspace_2.0-0
## [13] rstudioapi_0.13 stats4_4.0.3 robustbase_0.93-7
## [16] bayesplot_1.8.0 listenv_0.8.0 emmeans_1.5.3
## [19] RLRsim_3.1-6 coda_0.19-4 parallelly_1.23.0
## [22] vctrs_0.3.6 generics_0.1.0 TH.data_1.0-10
## [25] xfun_0.20 R6_2.5.0 markdown_1.1
## [28] gamm4_0.2-6 projpred_2.0.2 assertthat_0.2.1
## [31] promises_1.1.1 scales_1.1.1 multcomp_1.4-15
## [34] nnet_7.3-14 gtable_0.3.0 globals_0.14.0
## [37] processx_3.4.5 sandwich_3.0-0 rlang_0.4.10
## [40] splines_4.0.3 broom_0.7.3 inline_0.3.17
## [43] yaml_2.2.1 reshape2_1.4.4 abind_1.4-5
## [46] modelr_0.1.8 threejs_0.3.3 crosstalk_1.1.1
## [49] backports_1.2.1 httpuv_1.5.5 rsconnect_0.8.16
## [52] tools_4.0.3 ellipsis_0.3.1 ggridges_0.5.3
## [55] plyr_1.8.6 base64enc_0.1-3 ps_1.5.0
## [58] prettyunits_1.1.1 zoo_1.8-8 haven_2.3.1
## [61] fs_1.5.0 survey_4.0 magrittr_2.0.1
## [64] data.table_1.13.6 openxlsx_4.2.3 colourpicker_1.1.0
## [67] reprex_0.3.0 mvtnorm_1.1-1 sjmisc_2.8.6
## [70] matrixStats_0.57.0 hms_1.0.0 shinyjs_2.0.0
## [73] mime_0.9 evaluate_0.14 xtable_1.8-4
## [76] shinystan_2.5.0 pbkrtest_0.5-0.1 rio_0.5.16
## [79] sjstats_0.18.1 readxl_1.3.1 ggeffects_1.0.1
## [82] rstantools_2.1.1 compiler_4.0.3 V8_3.4.0
## [85] crayon_1.3.4 minqa_1.2.4 htmltools_0.5.1
## [88] mgcv_1.8-33 later_1.1.0.1 RcppParallel_5.0.2
## [91] lubridate_1.7.9.2 DBI_1.1.0 sjlabelled_1.1.7
## [94] dbplyr_1.4.2 MASS_7.3-53 boot_1.3-25
## [97] cli_2.2.0 mitools_2.4 parallel_4.0.3
## [100] insight_0.11.1 igraph_1.2.6 pkgconfig_2.0.3
## [103] numDeriv_2016.8-1.1 foreign_0.8-81 binom_1.1-1
## [106] xml2_1.3.2 dygraphs_1.1.1.6 estimability_1.3
## [109] rvest_0.3.6 callr_3.5.1 digest_0.6.27
## [112] parameters_0.10.1 fastGHQuad_1.0 rmarkdown_2.6
## [115] cellranger_1.1.0 curl_4.3 shiny_1.5.0
## [118] gtools_3.8.2 nloptr_1.2.2.2 lifecycle_0.2.0
## [121] nlme_3.1-151 jsonlite_1.7.2 fansi_0.4.1
## [124] pillar_1.4.7 loo_2.4.1 fastmap_1.0.1
## [127] httr_1.4.2 plotrix_3.7-8 DEoptimR_1.0-8
## [130] pkgbuild_1.2.0 survival_3.2-7 glue_1.4.2
## [133] xts_0.12.1 bayestestR_0.8.0 zip_2.1.1
## [136] shinythemes_1.1.2 iterators_1.0.13 stringi_1.5.3
## [139] performance_0.6.1# Read data ======================================
col_types <- cols(
labid = col_factor(),
subid = col_factor(),
subid_unique = col_factor(),
CDI.form = col_factor(),
CDI.nwords = col_integer(),
CDI.prop = col_number(),
CDI.agerange = col_factor(),
CDI.agedays = col_integer(),
CDI.agemin = col_integer(),
CDI.agemax = col_integer(),
vocab_nwords = col_integer(),
standardized.score.CDI = col_character(),
standardized.score.CDI.num = col_number(),
IDS_pref = col_number(),
language = col_factor(),
language_zone = col_factor(),
CDI.error = col_logical(),
Notes = col_character(),
trial_order = col_factor(),
method = col_factor(),
age_days = col_integer(),
age_mo = col_number(),
age_group = col_factor(),
nae = col_logical(),
gender = col_factor(),
second_session = col_logical()
)
data.total <- read_csv("data/02b_processed.csv", col_types = col_types)
# TODO: add saved resultsBefore moving on with the analysis, we have to ready the data by (a) checking for colinearity between z_age_months and CDI.z_age_months and correcting this if necessary, and (b) setting up the contrasts described in our data analysis.
1.1 Colinearity check
First, we run a Kappa test on the possibility of colinearity between z_age_months and CDI.z_age_months.
# Run kappa test
k.age_months <- model.matrix(~ z_age_months + CDI.z_age_months, data = data.total) %>%
kappa(exact = T)With a value of 3.1953332, we do not have a colinearity issue and can proceed with the data analysis as planned (The criteria of indicating colinearity is that kappa > 10).
1.2 Contrast Setups
We need gender as an effect-coded factor, and method as a deviation-coded factor. This is achieved in R by using the contr.sum() function with the number of levels for each factor. Notably, when subsetting the UK sample, only two levels of method out of the three in total were left.
# Set contrasts on the total dataset =============
contrasts(data.total$gender) <- contr.sum(2)
contrasts(data.total$method) <- contr.sum(3)
# Create sub-datasets, with contrasts ============
## NAE
data.nae <- data.total %>% subset(language_zone == "NAE") %>% droplevels()
contrasts(data.nae$gender) <- contr.sum(2)
contrasts(data.nae$method) <- contr.sum(3)
## UK
data.uk <- data.total %>% subset(language_zone == "British") %>% droplevels()
contrasts(data.uk$gender) <- contr.sum(2)
contrasts(data.uk$method) <- contr.sum(2) #note that UK sample has only 2 levels, so sum of zero contrasts set to 2 levels
## Other
data.other <- data.total %>% subset(language_zone == "Other") %>% droplevels()
contrasts(data.other$gender) <- contr.sum(2)
contrasts(data.other$method) <- contr.sum(3)2 Descriptive Statistics
We first assess the amount of data we have overall per condition and their shape overall.
data.total %>%
group_by(language_zone, CDI.agerange, method, gender) %>%
summarise(N = n(), age = mean(CDI.agedays), sd = sd(CDI.agedays)) %>%
kable()## `summarise()` regrouping output by 'language_zone', 'CDI.agerange', 'method' (override with `.groups` argument)| language_zone | CDI.agerange | method | gender | N | age | sd |
|---|---|---|---|---|---|---|
| British | 18 | singlescreen | F | 24 | 551.2083 | 10.741950 |
| British | 18 | singlescreen | M | 24 | 550.5000 | 13.497182 |
| British | 18 | eyetracking | F | 9 | 548.5556 | 9.139353 |
| British | 18 | eyetracking | M | 11 | 547.8182 | 10.147100 |
| British | 24 | singlescreen | F | 18 | 741.6111 | 13.128948 |
| British | 24 | singlescreen | M | 15 | 745.0667 | 15.549307 |
| British | 24 | eyetracking | F | 13 | 731.5385 | 12.162890 |
| British | 24 | eyetracking | M | 5 | 737.8000 | 8.228001 |
| Other | 18 | singlescreen | F | 11 | 541.5455 | 7.160498 |
| Other | 18 | singlescreen | M | 13 | 544.3077 | 6.663140 |
| Other | 18 | eyetracking | F | 26 | 556.0769 | 14.133430 |
| Other | 18 | eyetracking | M | 30 | 560.8333 | 16.128970 |
| Other | 18 | hpp | F | 38 | 549.0526 | 10.812774 |
| Other | 18 | hpp | M | 38 | 555.1579 | 13.664961 |
| Other | 24 | singlescreen | F | 10 | 731.3000 | 14.802777 |
| Other | 24 | singlescreen | M | 12 | 721.1667 | 13.888081 |
| Other | 24 | eyetracking | F | 28 | 730.1786 | 9.996494 |
| Other | 24 | eyetracking | M | 25 | 730.1600 | 7.711896 |
| Other | 24 | hpp | F | 45 | 730.7333 | 17.357733 |
| Other | 24 | hpp | M | 35 | 730.5143 | 15.849237 |
| NAE | 18 | singlescreen | F | 19 | 554.9474 | 20.780530 |
| NAE | 18 | singlescreen | M | 14 | 581.2143 | 24.925052 |
| NAE | 18 | eyetracking | F | 1 | 549.0000 | NA |
| NAE | 18 | hpp | F | 56 | 560.9821 | 21.584920 |
| NAE | 18 | hpp | M | 57 | 559.6140 | 21.163272 |
| NAE | 24 | singlescreen | F | 23 | 737.1739 | 26.604258 |
| NAE | 24 | singlescreen | M | 20 | 739.6000 | 21.352678 |
| NAE | 24 | eyetracking | F | 2 | 756.5000 | 34.648232 |
| NAE | 24 | eyetracking | M | 1 | 745.0000 | NA |
| NAE | 24 | hpp | F | 54 | 733.9630 | 27.476895 |
| NAE | 24 | hpp | M | 60 | 727.9500 | 27.409899 |
More detailed information about Descriptive Statistics
#number of lab
data.total %>%
select(labid, language_zone) %>%
unique() %>%
group_by(language_zone) %>%
count()## # A tibble: 3 x 2
## # Groups: language_zone [3]
## language_zone n
## <fct> <int>
## 1 British 4
## 2 Other 8
## 3 NAE 8data.total %>%
group_by(language_zone, CDI.agerange) %>%
summarize(N = n())## `summarise()` regrouping output by 'language_zone' (override with `.groups` argument)## # A tibble: 6 x 3
## # Groups: language_zone [3]
## language_zone CDI.agerange N
## <fct> <fct> <int>
## 1 British 18 68
## 2 British 24 51
## 3 Other 18 156
## 4 Other 24 155
## 5 NAE 18 147
## 6 NAE 24 160# age range in each age group and language zone
data.total %>%
select(subid, language_zone, CDI.agedays, CDI.agerange) %>%
unique() %>%
group_by(language_zone, CDI.agerange) %>%
summarize(age_min = (min(CDI.agedays)/365.25*12),
age_max = (max(CDI.agedays)/365.25*12))## `summarise()` regrouping output by 'language_zone' (override with `.groups` argument)## # A tibble: 6 x 4
## # Groups: language_zone [3]
## language_zone CDI.agerange age_min age_max
## <fct> <fct> <dbl> <dbl>
## 1 British 18 17.4 19.2
## 2 British 24 23.4 25.1
## 3 Other 18 17.4 19.5
## 4 Other 24 22.5 25.2
## 5 NAE 18 16.6 20.0
## 6 NAE 24 22.2 25.9We then assess the data per lab in terms of sample size and CDI score (vocabulary size, for consistency between language zones).
by_lab <- data.total %>%
group_by(labid, language_zone, CDI.agerange) %>%
mutate(tested = n_distinct(subid_unique)) %>%
select(labid, language_zone, CDI.agerange, tested, vocab_nwords) %>%
nest(scores = vocab_nwords) %>%
mutate(model = map(scores, ~ lm(vocab_nwords ~ 1, data = .x)),
ci = map(model, confint)) %>%
transmute(tested = tested,
mean = map_dbl(model, ~ coefficients(.x)[[1]]),
ci_lower = map_dbl(ci, 1),
ci_upper = map_dbl(ci, 2)) %>%
arrange(language_zone) %>%
rownames_to_column()
# TODO: find a way to group by language zone?
ggplot(by_lab,
aes(x = labid, colour = language_zone,
y = mean, ymin = ci_lower, ymax = ci_upper)) +
geom_linerange() +
geom_point(aes(size = tested)) +
facet_grid(cols = vars(CDI.agerange), scales = "free") + coord_flip(ylim = c(0, 500)) +
xlab("Lab") + ylab("Vocabulary size") +
scale_colour_brewer(palette = "Dark2", name = "Language zone",
breaks = c("British", "NAE", "Other")) +
scale_size_continuous(name = "N", breaks = function(x) c(min(x), mean(x), max(x))) +
theme(legend.position = "bottom")
3 Sample Theory Based Statistics
3.1 Simple Correlation
First, we want to assess quickly if there is a direct correlation between IDS preference and CDI score, computing a Pearson#’s product-moment correlation. We use standardized CDI scores for the North American sample, and raw scores for the British sample. [NOTE FROM MELANIE: I’m not sure this makes sense to do with the raw UK scores - since we’re collapsing across 18 and 24 month data and the data aren#’t standardized by age.].
# Statistics =====================================
## North American Sample
test.pearson.nae <- cor.test(data.nae$IDS_pref,
data.nae$z_standardized_CDI,
alternative = "two.sided", method = "pearson")
test.pearson.nae##
## Pearson's product-moment correlation
##
## data: data.nae$IDS_pref and data.nae$z_standardized_CDI
## t = -0.91847, df = 305, p-value = 0.3591
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## -0.16349833 0.05977293
## sample estimates:
## cor
## -0.05251901## UK Sample
test.pearson.uk <- cor.test(data.uk$IDS_pref,
data.uk$z_vocab_nwords,
alternative = "two.sided", method = "pearson")
test.pearson.uk##
## Pearson's product-moment correlation
##
## data: data.uk$IDS_pref and data.uk$z_vocab_nwords
## t = 1.0327, df = 117, p-value = 0.3039
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## -0.08643398 0.27040987
## sample estimates:
## cor
## 0.09504018Plots for correlation
## North American Sample
### Get correlation value for annotation
cor_text <- "paste(italic(R)^2, \" =\")"
cor_value <- round(test.pearson.nae$estimate, 3)
### Build plot
plot.pearson.nae <- data.nae %>%
ggplot(aes(x = IDS_pref,
y = standardized.score.CDI.num)) +
xlab("IDS preference") + ylab("Standardized CDI score") +
geom_point() +
geom_smooth(method = lm) +
annotate("text", x = -.9, y = 51, parse = T, size = 4,
label = paste(cor_text, cor_value, sep = "~"))
## UK Sample
cor_value <- round(test.pearson.uk$estimate, 3)
plot.pearson.uk <- data.uk %>%
ggplot(aes(x = IDS_pref,
y = vocab_nwords)) +
xlab("IDS preference") + ylab("Vocabulary size (in words)") +
geom_point() +
geom_smooth(method = lm) +
annotate("text", x = .8, y = 150, parse = T, size = 4,
label = paste(cor_text, cor_value, sep = "~"))
# Global plot
plot.pearson <- ggarrange(plot.pearson.nae, plot.pearson.uk, ncol = 2)## `geom_smooth()` using formula 'y ~ x'
## `geom_smooth()` using formula 'y ~ x'
plot.pearson
ggsave("plots/corr_plot.pdf", plot.pearson,
units = "mm", width = 180, height = 100, dpi = 1000)We see no obvious direct link between IDS prefernce and CDI score here. However, an effect might appear once we take into account various factors that might interact with IDS preference and/or CDI score. We can also first enhance these plots with information about the age group at which infants were tested (18- or 24-month-old), using vocabulary size to better compare the NAE and UK samples.
plot.age_group <- data.total %>%
subset(language_zone != "Other") %>%
droplevels() %>%
ggplot(aes(x = IDS_pref,
y = vocab_nwords,
colour = CDI.agerange)) +
facet_wrap(vars(language_zone),
labeller = as_labeller(c("British" = "UK samples",
"NAE" = "North Amercian samples"))) +
xlab("Standardized IDS prefence score") + ylab("Vocabulary size (in words)") + theme(legend.position = "top") +
geom_point() +
geom_smooth(method = lm) +
scale_colour_brewer(palette = "Dark2", name = "Age group",
breaks = c("18", "24"),
labels = c("18mo", "24m"))
ggsave("plots/scatter_age.pdf", plot.age_group,
units = "mm", width = 180, height = 100, dpi = 1000)## `geom_smooth()` using formula 'y ~ x'(plot.age_group)## `geom_smooth()` using formula 'y ~ x'
3.2 Mixed-Effects Model by Language Zone
Here, we run a mixed-effects model including only theoretically motivated effects, as described in the pre-registration. We start with the full model bellow, simplifying the random effects structure until it converges.
3.2.1 NAE full model
# Run models =====================================
## NAE
data.nae$centered_IDS_pref <- scale(data.nae$IDS_pref, center = T, scale = F)
lmer.full.nae <- lmer(standardized.score.CDI.num ~ CDI.z_age_months + gender +
z_age_months + method + centered_IDS_pref +
centered_IDS_pref:method + centered_IDS_pref:CDI.z_age_months + centered_IDS_pref:z_age_months +
(1 | labid),
data = data.nae)
summary(lmer.full.nae)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula:
## standardized.score.CDI.num ~ CDI.z_age_months + gender + z_age_months +
## method + centered_IDS_pref + centered_IDS_pref:method + centered_IDS_pref:CDI.z_age_months +
## centered_IDS_pref:z_age_months + (1 | labid)
## Data: data.nae
##
## REML criterion at convergence: 2861
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.7649 -0.8396 -0.1049 0.7743 2.3582
##
## Random effects:
## Groups Name Variance Std.Dev.
## labid (Intercept) 78.28 8.847
## Residual 752.79 27.437
## Number of obs: 307, groups: labid, 8
##
## Fixed effects:
## Estimate Std. Error df t value
## (Intercept) 31.38428 6.13351 28.72420 5.117
## CDI.z_age_months 0.07662 0.55332 295.20197 0.138
## gender1 -1.07741 1.61552 291.53707 -0.667
## z_age_months 0.35211 0.56449 153.40385 0.624
## method1 8.62370 5.88516 96.93942 1.465
## method2 -16.48345 10.24023 225.60947 -1.610
## centered_IDS_pref -36.37414 23.08356 289.74626 -1.576
## method1:centered_IDS_pref 27.84822 23.61831 290.27716 1.179
## method2:centered_IDS_pref -58.54287 46.03821 289.28918 -1.272
## CDI.z_age_months:centered_IDS_pref 1.15753 1.67303 290.78963 0.692
## z_age_months:centered_IDS_pref 1.78245 1.73910 292.19209 1.025
## Pr(>|t|)
## (Intercept) 1.88e-05 ***
## CDI.z_age_months 0.890
## gender1 0.505
## z_age_months 0.534
## method1 0.146
## method2 0.109
## centered_IDS_pref 0.116
## method1:centered_IDS_pref 0.239
## method2:centered_IDS_pref 0.205
## CDI.z_age_months:centered_IDS_pref 0.490
## z_age_months:centered_IDS_pref 0.306
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) CDI.z__ gendr1 z_g_mn methd1 methd2 c_IDS_ m1:_ID m2:_ID
## CDI.z_g_mnt -0.079
## gender1 -0.063 0.031
## z_age_mnths 0.031 -0.043 -0.070
## method1 -0.562 -0.003 0.019 0.140
## method2 0.800 -0.045 -0.057 -0.010 -0.777
## cntrd_IDS_p 0.202 0.019 0.027 -0.037 -0.203 0.238
## mthd1:_IDS_ -0.184 -0.049 -0.014 0.019 0.231 -0.244 -0.929
## mthd2:_IDS_ 0.195 0.040 0.047 -0.012 -0.222 0.244 0.970 -0.972
## CDI.__:_IDS -0.012 -0.021 -0.030 -0.018 -0.043 0.018 -0.075 0.038 -0.072
## z_g_m:_IDS_ -0.029 0.044 0.131 0.035 -0.049 0.030 0.037 -0.063 0.121
## CDI.__:
## CDI.z_g_mnt
## gender1
## z_age_mnths
## method1
## method2
## cntrd_IDS_p
## mthd1:_IDS_
## mthd2:_IDS_
## CDI.__:_IDS
## z_g_m:_IDS_ -0.097robust_lmer.full.nae <- robustlmm::rlmer(standardized.score.CDI.num ~ CDI.z_age_months + gender +
z_age_months + method + centered_IDS_pref +
centered_IDS_pref:method + centered_IDS_pref:CDI.z_age_months + centered_IDS_pref:z_age_months +
(1 | labid),
data = data.nae)
summary(lmer.full.nae)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula:
## standardized.score.CDI.num ~ CDI.z_age_months + gender + z_age_months +
## method + centered_IDS_pref + centered_IDS_pref:method + centered_IDS_pref:CDI.z_age_months +
## centered_IDS_pref:z_age_months + (1 | labid)
## Data: data.nae
##
## REML criterion at convergence: 2861
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.7649 -0.8396 -0.1049 0.7743 2.3582
##
## Random effects:
## Groups Name Variance Std.Dev.
## labid (Intercept) 78.28 8.847
## Residual 752.79 27.437
## Number of obs: 307, groups: labid, 8
##
## Fixed effects:
## Estimate Std. Error df t value
## (Intercept) 31.38428 6.13351 28.72420 5.117
## CDI.z_age_months 0.07662 0.55332 295.20197 0.138
## gender1 -1.07741 1.61552 291.53707 -0.667
## z_age_months 0.35211 0.56449 153.40385 0.624
## method1 8.62370 5.88516 96.93942 1.465
## method2 -16.48345 10.24023 225.60947 -1.610
## centered_IDS_pref -36.37414 23.08356 289.74626 -1.576
## method1:centered_IDS_pref 27.84822 23.61831 290.27716 1.179
## method2:centered_IDS_pref -58.54287 46.03821 289.28918 -1.272
## CDI.z_age_months:centered_IDS_pref 1.15753 1.67303 290.78963 0.692
## z_age_months:centered_IDS_pref 1.78245 1.73910 292.19209 1.025
## Pr(>|t|)
## (Intercept) 1.88e-05 ***
## CDI.z_age_months 0.890
## gender1 0.505
## z_age_months 0.534
## method1 0.146
## method2 0.109
## centered_IDS_pref 0.116
## method1:centered_IDS_pref 0.239
## method2:centered_IDS_pref 0.205
## CDI.z_age_months:centered_IDS_pref 0.490
## z_age_months:centered_IDS_pref 0.306
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) CDI.z__ gendr1 z_g_mn methd1 methd2 c_IDS_ m1:_ID m2:_ID
## CDI.z_g_mnt -0.079
## gender1 -0.063 0.031
## z_age_mnths 0.031 -0.043 -0.070
## method1 -0.562 -0.003 0.019 0.140
## method2 0.800 -0.045 -0.057 -0.010 -0.777
## cntrd_IDS_p 0.202 0.019 0.027 -0.037 -0.203 0.238
## mthd1:_IDS_ -0.184 -0.049 -0.014 0.019 0.231 -0.244 -0.929
## mthd2:_IDS_ 0.195 0.040 0.047 -0.012 -0.222 0.244 0.970 -0.972
## CDI.__:_IDS -0.012 -0.021 -0.030 -0.018 -0.043 0.018 -0.075 0.038 -0.072
## z_g_m:_IDS_ -0.029 0.044 0.131 0.035 -0.049 0.030 0.037 -0.063 0.121
## CDI.__:
## CDI.z_g_mnt
## gender1
## z_age_mnths
## method1
## method2
## cntrd_IDS_p
## mthd1:_IDS_
## mthd2:_IDS_
## CDI.__:_IDS
## z_g_m:_IDS_ -0.097full.nae_pvalue <- anova(lmer.full.nae) %>%
as_tibble(rownames = "Parameter") #this gives us the Type III p values
# ==========
# Sequentially removed random effects:
# IDS_pref:z_age_months
# IDS_pref:CDI.z_age_months
# IDS_pref:method
# IDS_pref
# method
# z_age_months
# gender
# CDI.z_age_months
#==========3.2.1.1 (Optional) Checking mixed-model assumptions. We will check the following:
- Linearlity
- Normality of the residuals
- Homoscedasticity of residuals
- No autocorrelation
- No multicollinearity
#First, check linearity
# data.nae$resid <- residuals(lmer.full.nae)
#
# plot(data.nae$resid, data.nae$standardized.score.CDI)
#Second, check normality
plot_model(lmer.full.nae, type = 'diag') ## we do have right-skewed normality of residuals## Warning in checkMatrixPackageVersion(): Package version inconsistency detected.
## TMB was built with Matrix version 1.2.18
## Current Matrix version is 1.3.2
## Please re-install 'TMB' from source using install.packages('TMB', type = 'source') or ask CRAN for a binary version of 'TMB' matching CRAN's 'Matrix' package## [[1]]## `geom_smooth()` using formula 'y ~ x'
##
## [[2]]
## [[2]]$labid## `geom_smooth()` using formula 'y ~ x'
##
##
## [[3]]
##
## [[4]]## `geom_smooth()` using formula 'y ~ x'
#Third, check autocorrelation
re_run_lme.full.nae <- nlme::lme(standardized.score.CDI.num ~ CDI.z_age_months + gender +
z_age_months + method + centered_IDS_pref +
centered_IDS_pref:method + centered_IDS_pref:CDI.z_age_months + centered_IDS_pref:z_age_months, random = ~1 | labid,
method = "REML",
data = data.nae, na.action = na.exclude)
plot(nlme::ACF(re_run_lme.full.nae, resType = "normalized")) #there is no sign for autocorrelation
#Lastly, check multi-collinearity
car::vif(lmer.full.nae) #we do see a multicollineartiy for the IDS preference variable, even though we have centered the IDS preference score. It is probably related to the number the participating labs (as this is the group level that we are controlling) and how we entered interaction between IDS preference and other variables (that lack variability in the current sample). We need to keep IDS preference in the model as exploring the relationship between IDS preference and CDI score is the key research question in the paper.## GVIF Df GVIF^(1/(2*Df))
## CDI.z_age_months 1.018050 1 1.008985
## gender 1.052996 1 1.026156
## z_age_months 1.076038 1 1.037322
## method 1.147009 2 1.034884
## centered_IDS_pref 22.501425 1 4.743567
## method:centered_IDS_pref 25.126190 2 2.238884
## CDI.z_age_months:centered_IDS_pref 1.032456 1 1.016099
## z_age_months:centered_IDS_pref 1.306986 1 1.1432353.2.2 UK full model
lmer.full.uk <- lmer(vocab_nwords ~ CDI.z_age_months + gender +
z_age_months + method + IDS_pref +
IDS_pref:method + IDS_pref:CDI.z_age_months + IDS_pref:z_age_months +
(1 | labid),
data = data.uk)
summary(lmer.full.uk)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: vocab_nwords ~ CDI.z_age_months + gender + z_age_months + method +
## IDS_pref + IDS_pref:method + IDS_pref:CDI.z_age_months +
## IDS_pref:z_age_months + (1 | labid)
## Data: data.uk
##
## REML criterion at convergence: 1342.7
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.0447 -0.6691 -0.0987 0.5147 3.5513
##
## Random effects:
## Groups Name Variance Std.Dev.
## labid (Intercept) 148.7 12.20
## Residual 7856.7 88.64
## Number of obs: 119, groups: labid, 4
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 175.103 11.540 3.162 15.174 0.000464 ***
## CDI.z_age_months 28.628 2.824 109.795 10.138 < 2e-16 ***
## gender1 17.578 8.623 109.015 2.039 0.043906 *
## z_age_months -6.034 3.042 101.094 -1.984 0.050013 .
## method1 -16.483 12.871 4.413 -1.281 0.263499
## IDS_pref 4.201 23.986 109.837 0.175 0.861288
## method1:IDS_pref 4.413 27.446 109.729 0.161 0.872550
## CDI.z_age_months:IDS_pref 0.600 7.877 108.111 0.076 0.939422
## z_age_months:IDS_pref 1.143 8.291 109.890 0.138 0.890629
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) CDI.z__ gendr1 z_g_mn methd1 IDS_pr m1:IDS CDI.__:
## CDI.z_g_mnt 0.148
## gender1 -0.023 -0.099
## z_age_mnths -0.005 -0.106 -0.001
## method1 -0.287 -0.075 -0.029 0.454
## IDS_pref -0.314 -0.092 -0.118 0.059 0.179
## mthd1:IDS_p 0.163 0.122 -0.012 -0.119 -0.264 -0.157
## CDI.__:IDS_ -0.126 -0.328 0.076 0.063 0.103 0.112 -0.195
## z_g_mn:IDS_ 0.014 0.095 -0.232 -0.360 -0.106 -0.096 0.493 -0.280full.uk_pvalue <- anova(lmer.full.uk) %>%
as_tibble(rownames = "Parameter") #this gives us the Type III p values
#==========
# Sequentially removed random effects:
# IDS_pref:z_age_months
# IDS_pref:CDI.z_age_months
# IDS_pref:method
# IDS_pref
# method
# z_age_months
# gender
# CDI.z_age_months3.2.2.1 (Optional) Checking mixed-model assumptions. We will check the following:
- Linearlity
- Normality of the residuals
- Homoscedasticity of residuals
- No autocorrelation
- No multicollinearity
#First, check linearity. The plot looks linear
data.uk$resid <- residuals(lmer.full.uk)
plot(data.uk$resid, data.uk$vocab_nwords)
#Second, check normality
plot_model(lmer.full.uk, type = 'diag') ## we do have right-skewed normality of residuals## [[1]]## `geom_smooth()` using formula 'y ~ x'
##
## [[2]]
## [[2]]$labid## `geom_smooth()` using formula 'y ~ x'
##
##
## [[3]]
##
## [[4]]## `geom_smooth()` using formula 'y ~ x'
#Third, check autocorrelation
re_run_lme.full.uk <- nlme::lme(vocab_nwords ~ CDI.z_age_months + gender +
z_age_months + method + IDS_pref +
IDS_pref:method + IDS_pref:CDI.z_age_months + IDS_pref:z_age_months, random = ~1 | labid,
method = "REML",
data = data.nae, na.action = na.exclude)
plot(nlme::ACF(re_run_lme.full.uk, resType = "normalized")) #there is no sign for autocorrelation
#Lastly, check multi-collinearity
car::vif(lmer.full.uk) #no problem for multicollinearlity## CDI.z_age_months gender z_age_months
## 1.155355 1.118400 1.553116
## method IDS_pref method:IDS_pref
## 1.445581 1.078478 1.515016
## CDI.z_age_months:IDS_pref z_age_months:IDS_pref
## 1.227720 1.819286We now want to check the statistical power of significant effects, and discard any models with significant effects that do not reach 80% power. This however leads to too many warnings of singularity issues on the model updates inherent to the simr power simulations, hence we cannot obtain satisfactory power estimates as pre-registered.
AST: Note that we don’t have any IV(s) that turned out to be significant in the Full NAE model. So we won’t run the power analysis check. For the UK full model, there are two statistically significant IV: CDI_age and gender. The post hoc power check suggested that we have high power in detecting the effect of CDI_age but not gender. Note that gender has a smaller effect size to begin with, so this may partially explain why we have less power in detecting it in the model. As there can be a number of different factors that determines the posthoc power, we decided not to remove gender in the model based on posthoc power analysis check.
check_pwr_uk_cdi_age <- simr::powerSim(lmer.full.uk, test = fixed("CDI.z_age_months", method = "z"), seed = 2, nsim = 1000, alpha = 0.05) #specify that Gender is the ixed effect that we are looking into
check_pwr_uk_cdi_age
check_pwr_uk_gender <- simr::powerSim(lmer.full.uk, test = fixed("gender", method = "z"), seed = 3, nsim = 1000, alpha = 0.05) #specify that Gender is the ixed effect that we are looking into
check_pwr_uk_gender3.2.3 Combined Sample
For this combined analysis, we first need to restrain the age range for the NAE sample (previously ±2 months, now ±1.5 months).
# Create dataset with British and NAE only
data.uk_nae <- data.total %>%
subset(language_zone %in% c("British", "NAE")) %>%
mutate(CDI.agemin = ifelse(language_zone == "NAE",
CDI.agemin + round(.5*365.25/12),
CDI.agemin),
CDI.agemax = ifelse(language_zone == "NAE",
CDI.agemax - round(.5*365.25/12),
CDI.agemax)) %>%
subset(!(CDI.agedays < CDI.agemin | CDI.agedays > CDI.agemax)) %>%
droplevels()
# Create contrasts for analysis
contrasts(data.uk_nae$gender) <- contr.sum(2)
contrasts(data.uk_nae$method) <- contr.sum(3)
contrasts(data.uk_nae$language_zone) <- contr.sum(2)We can then run the planned combined analysis adding the main effect and interactions of language_zone.
lmer.full.uk_nae <- lmer(CDI.prop ~ CDI.z_age_months + language_zone + gender +
z_age_months + method + IDS_pref + IDS_pref:language_zone +
IDS_pref:method + IDS_pref:CDI.z_age_months + IDS_pref:z_age_months +
(1 + CDI.z_age_months | labid),
data = data.uk_nae)
summary(lmer.full.uk_nae)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: CDI.prop ~ CDI.z_age_months + language_zone + gender + z_age_months +
## method + IDS_pref + IDS_pref:language_zone + IDS_pref:method +
## IDS_pref:CDI.z_age_months + IDS_pref:z_age_months + (1 +
## CDI.z_age_months | labid)
## Data: data.uk_nae
##
## REML criterion at convergence: -16.4
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.2182 -0.6666 -0.1781 0.5001 3.4898
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## labid (Intercept) 0.0027125 0.05208
## CDI.z_age_months 0.0008073 0.02841 0.40
## Residual 0.0439990 0.20976
## Number of obs: 401, groups: labid, 12
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 0.323365 0.022302 6.812745 14.499 2.26e-06 ***
## CDI.z_age_months 0.042960 0.009300 12.434234 4.619 0.000539 ***
## language_zone1 0.080976 0.027480 7.098782 2.947 0.021166 *
## gender1 0.037004 0.010822 372.652269 3.419 0.000697 ***
## z_age_months -0.001733 0.003994 289.264152 -0.434 0.664725
## method1 -0.010715 0.027902 12.016228 -0.384 0.707678
## method2 -0.001692 0.042360 10.268417 -0.040 0.968902
## IDS_pref -0.005438 0.038612 382.246227 -0.141 0.888062
## language_zone1:IDS_pref 0.068171 0.053566 374.217018 1.273 0.203927
## method1:IDS_pref -0.018179 0.051697 385.237670 -0.352 0.725301
## method2:IDS_pref -0.049967 0.080989 382.643250 -0.617 0.537625
## CDI.z_age_months:IDS_pref 0.008436 0.011007 383.996827 0.766 0.443883
## z_age_months:IDS_pref -0.003172 0.011006 376.867698 -0.288 0.773312
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1##
## Correlation matrix not shown by default, as p = 13 > 12.
## Use print(x, correlation=TRUE) or
## vcov(x) if you need itcombined.full.uk_nae_pvalue <- anova(lmer.full.uk_nae) %>%
as_tibble(rownames = "Parameter") #this gives us the Type III p values
#==========
# Sequentially removed random effects:
# IDS_pref:z_age_months
# IDS_pref:CDI.z_age_months
# IDS_pref:method
# IDS_pref:language_zone
# IDS_pref
# method
# z_age_months
# gender
# language_zone
#==========3.2.3.1 (Optional) Checking mixed-model assumptions
- Linearlity
- Normality of the residuals
- Homoscedasticity of residuals
- No autocorrelation
- No multicollinearity
#First, check linearity. The plot looks linear
data.uk_nae$resid <- residuals(lmer.full.uk_nae)
plot(data.uk_nae$resid, data.uk_nae$CDI.prop)
#Second, check normality
plot_model(lmer.full.uk_nae, type = 'diag') ## we do have right-skewed normality of residuals## [[1]]## `geom_smooth()` using formula 'y ~ x'
##
## [[2]]
## [[2]]$labid## `geom_smooth()` using formula 'y ~ x'
##
##
## [[3]]
##
## [[4]]## `geom_smooth()` using formula 'y ~ x'
#Third, check autocorrelation
re_run_lme.full.uk_nae <- nlme::lme(CDI.prop ~ CDI.z_age_months + language_zone + gender +
z_age_months + method + IDS_pref + IDS_pref:language_zone +
IDS_pref:method + IDS_pref:CDI.z_age_months + IDS_pref:z_age_months,
random = ~CDI.z_age_months | labid,
method = "REML",
data = data.uk_nae, na.action = na.exclude)
plot(nlme::ACF(re_run_lme.full.uk_nae, resType = "normalized")) #there is no sign for autocorrelation
#Lastly, check multi-collinearity
car::vif(lmer.full.uk_nae) #no problem for multicollinearlity## GVIF Df GVIF^(1/(2*Df))
## CDI.z_age_months 1.032963 1 1.016348
## language_zone 1.979191 1 1.406837
## gender 1.039619 1 1.019617
## z_age_months 1.328619 1 1.152657
## method 2.447762 2 1.250813
## IDS_pref 1.449858 1 1.204101
## language_zone:IDS_pref 2.900147 1 1.702982
## method:IDS_pref 3.854846 2 1.401205
## CDI.z_age_months:IDS_pref 1.078314 1 1.038419
## z_age_months:IDS_pref 1.509419 1 1.228584We then compute \(p\)-values, but leave out power estimates for those \(p\)-values as above. Again, we have a lot of singular fit issues for the power checks and decided not to remove parameters based on posthoc power analysis.
check_pwr_combined_cdi_age <- simr::powerSim(lmer.full.uk_nae, test = fixed("CDI.z_age_months", method = "z"), seed = 2, nsim = 1000, alpha = 0.05) #specify that Gender is the ixed effect that we are looking into
check_pwr_combined_cdi_age
check_pwr_combined_lang_zone <- simr::powerSim(lmer.full.uk_nae, test = fixed("language_zone", method = "z"), seed = 3, nsim = 1000, alpha = 0.05) #specify that Gender is the ixed effect that we are looking into
check_pwr_combined_lang_zone
check_pwr_combined_gender <- simr::powerSim(lmer.full.uk_nae, test = fixed("gender", method = "z"), seed = 3, nsim = 1000, alpha = 0.05) #specify that Gender is the ixed effect that we are looking into
check_pwr_combined_gender4 Bayesian Statistics
Given that the `frequentist#’ tests above did not discover any expected significant effects, we now run a Bayesian analysis, specifically Bayes factors (or \(b\)-values) from model comparisons, to determine whether our data provide evidence for a true null effect or merely fail to provide evidence either way.
We first run models on the separate UK and NAE samples.
First of all, we need to set up the priors for the Bayeisan analysis. In the RR, we said that we would use a Cauchy distribution for priors over fixed effects with normally distributed priors over random effect parameters."
4.1 Set up priors
Set up a really weak and not too informative prior
prior <-c(set_prior("normal(10, 1)", class = "Intercept"), #we can't have negative CDI
set_prior("cauchy(0, 1)", class = "b"), #all IVs have the same prior..
set_prior("normal(0, 1)", class = "sd"))
prior_combined_mod <-c(set_prior("normal(10, 1)", class = "Intercept"),
set_prior("cauchy(0, 1)", class = "b"), #all IVs have the same prior..
set_prior("normal(0, 1)", class = "sd"),
set_prior("lkj(2)", class = "L")) #adding this as there is a random slope and we need a correlation matrix prior between slope and interceptSet up another prior that mainly changed the meaan for the fixed effect and see if there is big difference? Note: numbers changed but not to the extent that change BF conclusion.
prior_big <-c(set_prior("normal(10, 1)", class = "Intercept"), #we can't have negative CDI
set_prior("cauchy(10, 5)", class = "b"), #all IVs have the same prior..
set_prior("normal(0, 1)", class = "sd"))
prior_combined_big <-c(set_prior("normal(10, 1)", class = "Intercept"),
set_prior("cauchy(10, 5)", class = "b"), #all IVs have the same prior..
set_prior("normal(0, 1)", class = "sd"),
set_prior("lkj(2)", class = "L")) #adding this as there is a random slope and we need a correlation matrix prior between slope and interceptWe also try more informative prior.
4.2 NAE info prior
prior_nae_1 <-c(set_prior("normal(50, 1)", class = "Intercept"), #note that nae's DV is percentile score, but I will assume at 18 months, the average is at 50 percentile
set_prior("cauchy(5, 1)", class = "b", coef = "CDI.z_age_months"), #assume increase in 1 month may lead to 5 % percentile increase?
set_prior("cauchy(-15, 5)", class = "b", coef = "gender1"), # female is the baseline, normally females know more words than males
set_prior("cauchy(0, 1)", class = "b", coef = "IDS_pref"), #weak prior as don't know the effect of IDS preference to CDI score
set_prior("cauchy(0, 1)", class = "b", coef = "method1"), #weak prior as methods related to the MB1 study shouldn't influence CDI scores
set_prior("cauchy(0, 1)", class = "b", coef = "method2"), #weak prior as methods related to the MB1 study shouldn't influence CDI scores
set_prior("cauchy(0, 1)", class = "b", coef = "z_age_months"),
set_prior("normal(0, 1)", class = "sd"))
prior_nae_null_1 <-c(set_prior("normal(50, 1)", class = "Intercept"), #note that nae's DV is percentile score, but I will assume at 18 months, the average is at 50 percentile
set_prior("cauchy(5, 1)", class = "b", coef = "CDI.z_age_months"), #assume increase in 1 month may lead to 5 % percentile increase?
set_prior("cauchy(-15, 5)", class = "b", coef = "gender1"), # female is the baseline, normally females know more words than males
set_prior("cauchy(0, 1)", class = "b", coef = "method1"), #weak prior as methods related to the MB1 study shouldn't influence CDI scores
set_prior("cauchy(0, 1)", class = "b", coef = "method2"), #weak prior as methods related to the MB1 study shouldn't influence CDI scores
set_prior("cauchy(0, 1)", class = "b", coef = "z_age_months"),
set_prior("normal(0, 1)", class = "sd"))
prior_nae_2 <-c(set_prior("normal(50, 1)", class = "Intercept"), #note that nae's DV is percentile score, but I will assume at 18 months, the average is at 50 percentile
set_prior("cauchy(5, 1)", class = "b", coef = "CDI.z_age_months"), #assume increase in 1 month may lead to 5 % percentile increase?
set_prior("cauchy(-15, 5)", class = "b", coef = "gender1"), # female is the baseline, normally females know more words than males
set_prior("cauchy(0, 1)", class = "b", coef = "IDS_pref"), #weak prior as don't know the effect of IDS preference to CDI score
set_prior("cauchy(0, 1)", class = "b", coef = "method1"), #weak prior as methods related to the MB1 study shouldn't influence CDI scores
set_prior("cauchy(0, 1)", class = "b", coef = "method2"), #weak prior as methods related to the MB1 study shouldn't influence CDI scores
#note all interaction priors are weak as we don't really have much info to form a prior
set_prior("cauchy(0, 1)", class = "b", coef = "z_age_months"),
set_prior("cauchy(0, 1)", class = "b", coef = "CDI.z_age_months:IDS_pref"),
set_prior("cauchy(0, 1)", class = "b", coef = "method1:IDS_pref"),
set_prior("cauchy(0, 1)", class = "b", coef = "method2:IDS_pref"),
set_prior("cauchy(0, 1)", class = "b", coef = "z_age_months:IDS_pref"),
set_prior("normal(0, 1)", class = "sd"))
prior_nae_null_2 <-c(set_prior("normal(50, 1)", class = "Intercept"), #note that nae's DV is percentile score, but I will assume at 18 months, the average is at 50 percentile
set_prior("cauchy(5, 1)", class = "b", coef = "CDI.z_age_months"), #assume increase in 1 month may lead to 5 % percentile increase?
set_prior("cauchy(-15, 5)", class = "b", coef = "gender1"), # female is the baseline, normally females know more words than males
set_prior("cauchy(0, 1)", class = "b", coef = "IDS_pref"), #weak prior as don't know the effect of IDS preference to CDI score
set_prior("cauchy(0, 1)", class = "b", coef = "method1"), #weak prior as methods related to the MB1 study shouldn't influence CDI scores
set_prior("cauchy(0, 1)", class = "b", coef = "method2"), #weak prior as methods related to the MB1 study shouldn't influence CDI scores
#note all interaction priors are weak as we don't really have much info to form a prior
set_prior("cauchy(0, 1)", class = "b", coef = "z_age_months"),
set_prior("cauchy(0, 1)", class = "b", coef = "CDI.z_age_months:IDS_pref"),
set_prior("cauchy(0, 1)", class = "b", coef = "method1:IDS_pref"),
set_prior("cauchy(0, 1)", class = "b", coef = "method2:IDS_pref"),
set_prior("normal(0, 1)", class = "sd")) 4.3 UK info prior
prior_uk_1 <-c(set_prior("normal(40, 1)", class = "Intercept"), #note that UK's DV is raw score, based on wordbank (oxford CDI), the average (50 percentile) is around 41.6
set_prior("cauchy(15, 1)", class = "b", coef = "CDI.z_age_months"), #assume increase in 1 month may lead to 15 words increase (proxy from wordbank)
set_prior("cauchy(-10, 5)", class = "b", coef = "gender1"), # female is the baseline, normally females know more words than males. Wordbank oxford CDI suggested at 18 months (50 percentile), female knows approximately 11 words more than males
set_prior("cauchy(0, 1)", class = "b", coef = "IDS_pref"), #weak prior as don't know the effect of IDS preference to CDI score
set_prior("cauchy(0, 1)", class = "b", coef = "method1"), #weak prior as methods related to the MB1 study shouldn't influence CDI scores
set_prior("cauchy(0, 1)", class = "b", coef = "z_age_months"),
set_prior("normal(0, 1)", class = "sd"))
prior_uk_null_1 <-c(set_prior("normal(40, 1)", class = "Intercept"), #note that UK's DV is raw score, based on wordbank (oxford CDI), the average (50 percentile) is around 41.6
set_prior("cauchy(15, 1)", class = "b", coef = "CDI.z_age_months"), #assume increase in 1 month may lead to 15 words increase (proxy from wordbank)
set_prior("cauchy(-10, 5)", class = "b", coef = "gender1"), # female is the baseline, normally females know more words than males. Wordbank oxford CDI suggested at 18 months (50 percentile), female knows approximately 11 words more than males
set_prior("cauchy(0, 1)", class = "b", coef = "method1"), #weak prior as methods related to the MB1 study shouldn't influence CDI scores
set_prior("cauchy(0, 1)", class = "b", coef = "z_age_months"),
set_prior("normal(0, 1)", class = "sd"))
prior_uk_2 <-c(set_prior("normal(40, 1)", class = "Intercept"), #note that UK's DV is raw score, based on wordbank (oxford CDI), the average (50 percentile) is around 41.6
set_prior("cauchy(15, 1)", class = "b", coef = "CDI.z_age_months"), #assume increase in 1 month may lead to 15 words increase (proxy from wordbank)
set_prior("cauchy(-10, 5)", class = "b", coef = "gender1"), # female is the baseline, normally females know more words than males. Wordbank oxford CDI suggested at 18 months (50 percentile), female knows approximately 11 words more than males
set_prior("cauchy(0, 1)", class = "b", coef = "IDS_pref"), #weak prior as don't know the effect of IDS preference to CDI score
set_prior("cauchy(0, 1)", class = "b", coef = "method1"), #weak prior as methods related to the MB1 study shouldn't influence CDI scores
set_prior("cauchy(0, 1)", class = "b", coef = "z_age_months"),
#note all interaction priors are weak as we don't really have much info to form a prior
set_prior("cauchy(0, 1)", class = "b", coef = "CDI.z_age_months:IDS_pref"),
set_prior("cauchy(0, 1)", class = "b", coef = "method1:IDS_pref"),
set_prior("cauchy(0, 1)", class = "b", coef = "z_age_months:IDS_pref"),
set_prior("normal(0, 1)", class = "sd"))
prior_uk_null_2 <-c(set_prior("normal(40, 1)", class = "Intercept"), #note that UK's DV is raw score, based on wordbank (oxford CDI), the average (50 percentile) is around 41.6
set_prior("cauchy(15, 1)", class = "b", coef = "CDI.z_age_months"), #assume increase in 1 month may lead to 15 words increase (proxy from wordbank)
set_prior("cauchy(-10, 5)", class = "b", coef = "gender1"), # female is the baseline, normally females know more words than males. Wordbank oxford CDI suggested at 18 months (50 percentile), female knows approximately 11 words more than males
set_prior("cauchy(0, 1)", class = "b", coef = "IDS_pref"), #weak prior as don't know the effect of IDS preference to CDI score
set_prior("cauchy(0, 1)", class = "b", coef = "method1"), #weak prior as methods related to the MB1 study shouldn't influence CDI scores
set_prior("cauchy(0, 1)", class = "b", coef = "z_age_months"),
#note all interaction priors are weak as we don't really have much info to form a prior
set_prior("cauchy(0, 1)", class = "b", coef = "CDI.z_age_months:IDS_pref"),
set_prior("cauchy(0, 1)", class = "b", coef = "method1:IDS_pref"),
set_prior("normal(0, 1)", class = "sd")) 4.4 Combined info prior
# Note that NAE CDI WS total score is 680 and UK CDI total score is 416.
# NAE average scores for 50th percentile is 76.3, so it is 76.3/680 ~ 11% and UK CDI average scores for 50th percentile is 41.6, so it is 41.6/416 ~ 10%
prior_combined_full_info <-c(set_prior("normal(10, 1)", class = "Intercept"), # given the above info, we expect on average, proportion score is 10%
set_prior("cauchy(3, 1)", class = "b", coef = "CDI.z_age_months"), #assume increase in 1 month may lead to 15 words increase in UK cdi and lead to 28 words increase in NAE (proxy from wordbank). In general, 15/416 = 3% and 28/680 = 4%
set_prior("cauchy(-3, 1)", class = "b", coef = "gender1"), # female is the baseline, normally females know more words than males. Wordbank oxford CDI suggested at 18 months (50 percentile), female knows approximately 11 words more than males, 11/416 ~ 2.6% and Wordbank NAE suggested females knows approximately 27 words than males, 27/680 ~ 3.9%
set_prior("cauchy(0, 1)", class = "b", coef = "language_zone1"), #weak prior as don't know the effect of language zone to CDI score
set_prior("cauchy(0, 1)", class = "b", coef = "IDS_pref"), #weak prior as don't know the effect of IDS preference to CDI score
set_prior("cauchy(0, 1)", class = "b", coef = "method1"), #weak prior as methods related to the MB1 study shouldn't influence CDI scores
set_prior("cauchy(0, 1)", class = "b", coef = "method2"), #weak prior as methods related to the MB1 study shouldn't influence CDI scores
set_prior("cauchy(0, 1)", class = "b", coef = "z_age_months"),
#note all interaction priors are weak as we don't really have much info to form a prior
set_prior("cauchy(0, 1)", class = "b", coef = "CDI.z_age_months:IDS_pref"),
set_prior("cauchy(0, 1)", class = "b", coef = "language_zone1:IDS_pref"),
set_prior("cauchy(0, 1)", class = "b", coef = "method1:IDS_pref"),
set_prior("cauchy(0, 1)", class = "b", coef = "method2:IDS_pref"),
set_prior("cauchy(0, 1)", class = "b", coef = "z_age_months:IDS_pref"),
set_prior("normal(0, 1)", class = "sd"),
set_prior("lkj(2)", class = "L")) #adding this as there is a random slope and we need a correlation matrix prior between slope and intercept
prior_combined_null_info <-c(set_prior("normal(10, 1)", class = "Intercept"), # given the above info, we expect on average, proportion score is 10%
set_prior("cauchy(3, 1)", class = "b", coef = "CDI.z_age_months"), #assume increase in 1 month may lead to 15 words increase in UK cdi and lead to 28 words increase in NAE (proxy from wordbank). In general, 15/416 = 3% and 28/680 = 4%
set_prior("cauchy(-3, 1)", class = "b", coef = "gender1"), # female is the baseline, normally females know more words than males. Wordbank oxford CDI suggested at 18 months (50 percentile), female knows approximately 11 words more than males, 11/416 ~ 2.6% and Wordbank NAE suggested females knows approximately 27 words than males, 27/680 ~ 3.9%
set_prior("cauchy(0, 1)", class = "b", coef = "language_zone1"), #weak prior as don't know the effect of language zone to CDI score
set_prior("cauchy(0, 1)", class = "b", coef = "IDS_pref"), #weak prior as don't know the effect of IDS preference to CDI score
set_prior("cauchy(0, 1)", class = "b", coef = "method1"), #weak prior as methods related to the MB1 study shouldn't influence CDI scores
set_prior("cauchy(0, 1)", class = "b", coef = "method2"), #weak prior as methods related to the MB1 study shouldn't influence CDI scores
set_prior("cauchy(0, 1)", class = "b", coef = "z_age_months"),
#note all interaction priors are weak as we don't really have much info to form a prior
set_prior("cauchy(0, 1)", class = "b", coef = "CDI.z_age_months:IDS_pref"),
set_prior("cauchy(0, 1)", class = "b", coef = "method1:IDS_pref"),
set_prior("cauchy(0, 1)", class = "b", coef = "method2:IDS_pref"),
set_prior("cauchy(0, 1)", class = "b", coef = "z_age_months:IDS_pref"),
set_prior("normal(0, 1)", class = "sd"),
set_prior("lkj(2)", class = "L")) #adding this as there is a random slope and we need a correlation matrix prior between slope and interceptIn the following analysis, we follow what we have proposed in the RR about calculating the Bayes factor for the three variables of interest: 1) The main effect IDS preference on CDI scores 2) The interaction effect between IDS preference and IDS test age on CDI scores 3) The interaction effect between language_zone (i.e., dialect in RR) and IDS preference on CDI scores
For the first two variables, we run it on the NAE and UK model separately. For the third variable, we run in on the combined UK & NAE model.
To calcualte Bayes factor for the main effect IDS preference on CDI scores, we need to two models: a full model that contains the main effect of variables of interest, and the other null model that does not contain the main effect of interest. As all the interaction terms in the models contains IDS preference, in the following, we only include main effects in the full and null models for comparisons.
4.5 NAE model (main effect of IDS preference)
Sys.setenv(R_FUTURE_RNG_ONMISUSE = "ignore") #this line of code is needed to get rid of warning message due to the "future" package. https://github.com/satijalab/seurat/issues/3622
lmer.bf_full_1.nae <- brm(standardized.score.CDI.num ~ CDI.z_age_months + gender +
z_age_months + method + IDS_pref +
(1 | labid),
data = data.nae,
family = gaussian,
prior = prior_nae_1,
iter = 4000,
warmup = 2000,
chains = 4,
save_pars = save_pars(all = T),
control = list(adapt_delta = .95),
seed=456)## Compiling Stan program...## Start sampling##
## SAMPLING FOR MODEL '3b2c8dc852b93f8f4a16392e0b6193a9' NOW (CHAIN 1).
## Chain 1:
## Chain 1: Gradient evaluation took 7e-05 seconds
## Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 0.7 seconds.
## Chain 1: Adjust your expectations accordingly!
## Chain 1:
## Chain 1:
## Chain 1: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 1: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 1: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 1: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 1: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 1: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 1: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 1: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 1: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 1: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 1: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 1: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 1:
## Chain 1: Elapsed Time: 0.742789 seconds (Warm-up)
## Chain 1: 0.769977 seconds (Sampling)
## Chain 1: 1.51277 seconds (Total)
## Chain 1:
##
## SAMPLING FOR MODEL '3b2c8dc852b93f8f4a16392e0b6193a9' NOW (CHAIN 2).
## Chain 2:
## Chain 2: Gradient evaluation took 4.9e-05 seconds
## Chain 2: 1000 transitions using 10 leapfrog steps per transition would take 0.49 seconds.
## Chain 2: Adjust your expectations accordingly!
## Chain 2:
## Chain 2:
## Chain 2: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 2: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 2: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 2: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 2: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 2: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 2: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 2: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 2: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 2: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 2: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 2: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 2:
## Chain 2: Elapsed Time: 0.741755 seconds (Warm-up)
## Chain 2: 0.472661 seconds (Sampling)
## Chain 2: 1.21442 seconds (Total)
## Chain 2:
##
## SAMPLING FOR MODEL '3b2c8dc852b93f8f4a16392e0b6193a9' NOW (CHAIN 3).
## Chain 3:
## Chain 3: Gradient evaluation took 5.4e-05 seconds
## Chain 3: 1000 transitions using 10 leapfrog steps per transition would take 0.54 seconds.
## Chain 3: Adjust your expectations accordingly!
## Chain 3:
## Chain 3:
## Chain 3: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 3: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 3: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 3: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 3: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 3: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 3: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 3: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 3: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 3: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 3: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 3: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 3:
## Chain 3: Elapsed Time: 0.69613 seconds (Warm-up)
## Chain 3: 0.749866 seconds (Sampling)
## Chain 3: 1.446 seconds (Total)
## Chain 3:
##
## SAMPLING FOR MODEL '3b2c8dc852b93f8f4a16392e0b6193a9' NOW (CHAIN 4).
## Chain 4:
## Chain 4: Gradient evaluation took 4.9e-05 seconds
## Chain 4: 1000 transitions using 10 leapfrog steps per transition would take 0.49 seconds.
## Chain 4: Adjust your expectations accordingly!
## Chain 4:
## Chain 4:
## Chain 4: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 4: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 4: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 4: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 4: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 4: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 4: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 4: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 4: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 4: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 4: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 4: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 4:
## Chain 4: Elapsed Time: 0.693584 seconds (Warm-up)
## Chain 4: 0.662766 seconds (Sampling)
## Chain 4: 1.35635 seconds (Total)
## Chain 4:summary(lmer.bf_full_1.nae)## Family: gaussian
## Links: mu = identity; sigma = identity
## Formula: standardized.score.CDI.num ~ CDI.z_age_months + gender + z_age_months + method + IDS_pref + (1 | labid)
## Data: data.nae (Number of observations: 307)
## Samples: 4 chains, each with iter = 4000; warmup = 2000; thin = 1;
## total post-warmup samples = 8000
##
## Group-Level Effects:
## ~labid (Number of levels: 8)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 2.18 1.07 0.17 4.19 1.00 2834 2875
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## Intercept 47.60 1.67 44.19 50.82 1.00 6735 4199
## CDI.z_age_months 0.15 0.57 -0.99 1.28 1.00 9418 5836
## gender1 -1.98 1.69 -5.38 1.24 1.00 9021 6193
## z_age_months -0.04 0.44 -0.88 0.82 1.00 8768 5974
## method1 -0.01 1.30 -2.70 2.72 1.00 6826 4206
## method2 -0.28 1.97 -4.88 3.47 1.00 6414 3069
## IDS_pref -0.57 2.04 -5.98 2.79 1.00 7243 2554
##
## Family Specific Parameters:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma 28.66 1.25 26.34 31.23 1.00 8001 5743
##
## Samples were drawn using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).lmer.bf_null_1.nae <- brm(standardized.score.CDI.num ~ CDI.z_age_months + gender +
z_age_months + method +
(1 | labid),
data = data.nae,
family = gaussian,
prior = prior_nae_null_1,
iter = 4000,
warmup = 2000,
chains = 4,
save_pars = save_pars(all = T),
control = list(adapt_delta = .95),
seed =123)## Compiling Stan program...
## Start sampling##
## SAMPLING FOR MODEL '32a715b1500064463a75f1c5a10e4ab8' NOW (CHAIN 1).
## Chain 1:
## Chain 1: Gradient evaluation took 5.9e-05 seconds
## Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 0.59 seconds.
## Chain 1: Adjust your expectations accordingly!
## Chain 1:
## Chain 1:
## Chain 1: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 1: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 1: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 1: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 1: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 1: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 1: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 1: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 1: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 1: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 1: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 1: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 1:
## Chain 1: Elapsed Time: 0.657162 seconds (Warm-up)
## Chain 1: 0.438973 seconds (Sampling)
## Chain 1: 1.09614 seconds (Total)
## Chain 1:
##
## SAMPLING FOR MODEL '32a715b1500064463a75f1c5a10e4ab8' NOW (CHAIN 2).
## Chain 2:
## Chain 2: Gradient evaluation took 9e-05 seconds
## Chain 2: 1000 transitions using 10 leapfrog steps per transition would take 0.9 seconds.
## Chain 2: Adjust your expectations accordingly!
## Chain 2:
## Chain 2:
## Chain 2: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 2: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 2: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 2: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 2: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 2: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 2: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 2: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 2: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 2: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 2: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 2: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 2:
## Chain 2: Elapsed Time: 0.71534 seconds (Warm-up)
## Chain 2: 0.801815 seconds (Sampling)
## Chain 2: 1.51716 seconds (Total)
## Chain 2:
##
## SAMPLING FOR MODEL '32a715b1500064463a75f1c5a10e4ab8' NOW (CHAIN 3).
## Chain 3:
## Chain 3: Gradient evaluation took 8.5e-05 seconds
## Chain 3: 1000 transitions using 10 leapfrog steps per transition would take 0.85 seconds.
## Chain 3: Adjust your expectations accordingly!
## Chain 3:
## Chain 3:
## Chain 3: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 3: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 3: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 3: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 3: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 3: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 3: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 3: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 3: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 3: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 3: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 3: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 3:
## Chain 3: Elapsed Time: 0.58087 seconds (Warm-up)
## Chain 3: 0.812227 seconds (Sampling)
## Chain 3: 1.3931 seconds (Total)
## Chain 3:
##
## SAMPLING FOR MODEL '32a715b1500064463a75f1c5a10e4ab8' NOW (CHAIN 4).
## Chain 4:
## Chain 4: Gradient evaluation took 9.3e-05 seconds
## Chain 4: 1000 transitions using 10 leapfrog steps per transition would take 0.93 seconds.
## Chain 4: Adjust your expectations accordingly!
## Chain 4:
## Chain 4:
## Chain 4: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 4: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 4: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 4: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 4: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 4: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 4: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 4: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 4: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 4: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 4: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 4: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 4:
## Chain 4: Elapsed Time: 0.706565 seconds (Warm-up)
## Chain 4: 0.890216 seconds (Sampling)
## Chain 4: 1.59678 seconds (Total)
## Chain 4:summary(lmer.bf_null_1.nae)## Family: gaussian
## Links: mu = identity; sigma = identity
## Formula: standardized.score.CDI.num ~ CDI.z_age_months + gender + z_age_months + method + (1 | labid)
## Data: data.nae (Number of observations: 307)
## Samples: 4 chains, each with iter = 4000; warmup = 2000; thin = 1;
## total post-warmup samples = 8000
##
## Group-Level Effects:
## ~labid (Number of levels: 8)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 2.18 1.07 0.17 4.14 1.00 3174 2530
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## Intercept 47.50 1.63 44.04 50.60 1.00 8004 4014
## CDI.z_age_months 0.13 0.57 -0.98 1.24 1.00 11260 5283
## gender1 -1.96 1.67 -5.24 1.30 1.00 12196 5268
## z_age_months -0.05 0.43 -0.90 0.80 1.00 11669 5684
## method1 0.00 1.30 -2.72 2.69 1.00 8248 4403
## method2 -0.30 1.94 -4.93 3.34 1.00 6939 3312
##
## Family Specific Parameters:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma 28.66 1.25 26.33 31.19 1.00 8247 5849
##
## Samples were drawn using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).Calculate the marginal log likelihoods for each hypothese using Bridge sampling
marloglik_full_1.nae <- bridge_sampler(lmer.bf_full_1.nae, silent = T)
marloglik_null_1.nae <- bridge_sampler(lmer.bf_null_1.nae, silent = T)Bayes factor
BF_full_1.nae <- bayes_factor(marloglik_full_1.nae, marloglik_null_1.nae)
BF_full_1.nae## Estimated Bayes factor in favor of x1 over x2: 0.914184.6 NAE model (interaction between IDS test age and IDS preference)
lmer.bf_full_2.nae <- brm(standardized.score.CDI.num ~ CDI.z_age_months + gender +
z_age_months + method + IDS_pref +
IDS_pref:method + IDS_pref:CDI.z_age_months + IDS_pref:z_age_months +
(1 | labid),
data = data.nae,
family = gaussian,
prior = prior_nae_2,
iter = 4000,
warmup = 2000,
chains = 4,
save_pars = save_pars(all = T),
control = list(adapt_delta = .95),
seed=890)## Compiling Stan program...## Start sampling##
## SAMPLING FOR MODEL 'c92e15777c97d2cdb32465f3d1c85d38' NOW (CHAIN 1).
## Chain 1:
## Chain 1: Gradient evaluation took 7.4e-05 seconds
## Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 0.74 seconds.
## Chain 1: Adjust your expectations accordingly!
## Chain 1:
## Chain 1:
## Chain 1: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 1: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 1: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 1: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 1: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 1: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 1: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 1: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 1: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 1: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 1: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 1: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 1:
## Chain 1: Elapsed Time: 0.987588 seconds (Warm-up)
## Chain 1: 0.771481 seconds (Sampling)
## Chain 1: 1.75907 seconds (Total)
## Chain 1:
##
## SAMPLING FOR MODEL 'c92e15777c97d2cdb32465f3d1c85d38' NOW (CHAIN 2).
## Chain 2:
## Chain 2: Gradient evaluation took 7.8e-05 seconds
## Chain 2: 1000 transitions using 10 leapfrog steps per transition would take 0.78 seconds.
## Chain 2: Adjust your expectations accordingly!
## Chain 2:
## Chain 2:
## Chain 2: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 2: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 2: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 2: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 2: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 2: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 2: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 2: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 2: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 2: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 2: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 2: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 2:
## Chain 2: Elapsed Time: 0.802676 seconds (Warm-up)
## Chain 2: 0.863793 seconds (Sampling)
## Chain 2: 1.66647 seconds (Total)
## Chain 2:
##
## SAMPLING FOR MODEL 'c92e15777c97d2cdb32465f3d1c85d38' NOW (CHAIN 3).
## Chain 3:
## Chain 3: Gradient evaluation took 6.9e-05 seconds
## Chain 3: 1000 transitions using 10 leapfrog steps per transition would take 0.69 seconds.
## Chain 3: Adjust your expectations accordingly!
## Chain 3:
## Chain 3:
## Chain 3: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 3: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 3: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 3: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 3: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 3: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 3: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 3: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 3: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 3: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 3: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 3: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 3:
## Chain 3: Elapsed Time: 0.812938 seconds (Warm-up)
## Chain 3: 0.726661 seconds (Sampling)
## Chain 3: 1.5396 seconds (Total)
## Chain 3:
##
## SAMPLING FOR MODEL 'c92e15777c97d2cdb32465f3d1c85d38' NOW (CHAIN 4).
## Chain 4:
## Chain 4: Gradient evaluation took 7.2e-05 seconds
## Chain 4: 1000 transitions using 10 leapfrog steps per transition would take 0.72 seconds.
## Chain 4: Adjust your expectations accordingly!
## Chain 4:
## Chain 4:
## Chain 4: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 4: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 4: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 4: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 4: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 4: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 4: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 4: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 4: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 4: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 4: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 4: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 4:
## Chain 4: Elapsed Time: 0.868436 seconds (Warm-up)
## Chain 4: 0.790037 seconds (Sampling)
## Chain 4: 1.65847 seconds (Total)
## Chain 4:summary(lmer.bf_full_2.nae)## Family: gaussian
## Links: mu = identity; sigma = identity
## Formula: standardized.score.CDI.num ~ CDI.z_age_months + gender + z_age_months + method + IDS_pref + IDS_pref:method + IDS_pref:CDI.z_age_months + IDS_pref:z_age_months + (1 | labid)
## Data: data.nae (Number of observations: 307)
## Samples: 4 chains, each with iter = 4000; warmup = 2000; thin = 1;
## total post-warmup samples = 8000
##
## Group-Level Effects:
## ~labid (Number of levels: 8)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 2.15 1.05 0.17 4.06 1.00 3214 3025
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept 47.44 1.74 43.88 50.66 1.00 7073
## CDI.z_age_months 0.08 0.60 -1.08 1.24 1.00 8998
## gender1 -1.90 1.66 -5.10 1.34 1.00 9048
## z_age_months -0.12 0.45 -1.00 0.77 1.00 8524
## method1 -0.01 1.35 -2.81 2.73 1.00 7236
## method2 -0.34 2.06 -5.14 3.45 1.00 6897
## IDS_pref -0.84 2.38 -7.36 2.66 1.00 6695
## method1:IDS_pref -0.00 2.00 -4.42 4.30 1.00 8656
## method2:IDS_pref 0.28 2.36 -4.27 6.11 1.00 7900
## CDI.z_age_months:IDS_pref 0.44 1.10 -1.55 2.94 1.00 8479
## z_age_months:IDS_pref 0.76 1.15 -1.16 3.45 1.00 7334
## Tail_ESS
## Intercept 4191
## CDI.z_age_months 5991
## gender1 5386
## z_age_months 6523
## method1 3651
## method2 3111
## IDS_pref 2868
## method1:IDS_pref 3268
## method2:IDS_pref 2693
## CDI.z_age_months:IDS_pref 4637
## z_age_months:IDS_pref 4296
##
## Family Specific Parameters:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma 28.66 1.26 26.32 31.19 1.00 7384 6053
##
## Samples were drawn using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).lmer.bf_null_2.nae <- brm(standardized.score.CDI.num ~ CDI.z_age_months + gender +
z_age_months + method + IDS_pref +
IDS_pref:method + IDS_pref:CDI.z_age_months +
(1 | labid),
data = data.nae,
family = gaussian,
prior = prior_nae_null_2,
iter = 4000,
warmup = 2000,
chains = 4,
save_pars = save_pars(all = T),
control = list(adapt_delta = .95),
seed =102)## Compiling Stan program...
## Start sampling##
## SAMPLING FOR MODEL '91b614181d9419722a0fc80a47971299' NOW (CHAIN 1).
## Chain 1:
## Chain 1: Gradient evaluation took 0.000262 seconds
## Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 2.62 seconds.
## Chain 1: Adjust your expectations accordingly!
## Chain 1:
## Chain 1:
## Chain 1: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 1: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 1: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 1: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 1: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 1: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 1: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 1: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 1: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 1: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 1: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 1: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 1:
## Chain 1: Elapsed Time: 0.83679 seconds (Warm-up)
## Chain 1: 0.982422 seconds (Sampling)
## Chain 1: 1.81921 seconds (Total)
## Chain 1:
##
## SAMPLING FOR MODEL '91b614181d9419722a0fc80a47971299' NOW (CHAIN 2).
## Chain 2:
## Chain 2: Gradient evaluation took 7.4e-05 seconds
## Chain 2: 1000 transitions using 10 leapfrog steps per transition would take 0.74 seconds.
## Chain 2: Adjust your expectations accordingly!
## Chain 2:
## Chain 2:
## Chain 2: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 2: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 2: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 2: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 2: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 2: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 2: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 2: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 2: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 2: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 2: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 2: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 2:
## Chain 2: Elapsed Time: 0.871272 seconds (Warm-up)
## Chain 2: 0.675621 seconds (Sampling)
## Chain 2: 1.54689 seconds (Total)
## Chain 2:
##
## SAMPLING FOR MODEL '91b614181d9419722a0fc80a47971299' NOW (CHAIN 3).
## Chain 3:
## Chain 3: Gradient evaluation took 7e-05 seconds
## Chain 3: 1000 transitions using 10 leapfrog steps per transition would take 0.7 seconds.
## Chain 3: Adjust your expectations accordingly!
## Chain 3:
## Chain 3:
## Chain 3: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 3: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 3: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 3: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 3: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 3: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 3: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 3: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 3: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 3: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 3: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 3: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 3:
## Chain 3: Elapsed Time: 0.962901 seconds (Warm-up)
## Chain 3: 0.55674 seconds (Sampling)
## Chain 3: 1.51964 seconds (Total)
## Chain 3:
##
## SAMPLING FOR MODEL '91b614181d9419722a0fc80a47971299' NOW (CHAIN 4).
## Chain 4:
## Chain 4: Gradient evaluation took 7.6e-05 seconds
## Chain 4: 1000 transitions using 10 leapfrog steps per transition would take 0.76 seconds.
## Chain 4: Adjust your expectations accordingly!
## Chain 4:
## Chain 4:
## Chain 4: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 4: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 4: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 4: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 4: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 4: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 4: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 4: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 4: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 4: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 4: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 4: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 4:
## Chain 4: Elapsed Time: 0.905299 seconds (Warm-up)
## Chain 4: 0.732245 seconds (Sampling)
## Chain 4: 1.63754 seconds (Total)
## Chain 4:summary(lmer.bf_null_2.nae)## Family: gaussian
## Links: mu = identity; sigma = identity
## Formula: standardized.score.CDI.num ~ CDI.z_age_months + gender + z_age_months + method + IDS_pref + IDS_pref:method + IDS_pref:CDI.z_age_months + (1 | labid)
## Data: data.nae (Number of observations: 307)
## Samples: 4 chains, each with iter = 4000; warmup = 2000; thin = 1;
## total post-warmup samples = 8000
##
## Group-Level Effects:
## ~labid (Number of levels: 8)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 2.19 1.06 0.18 4.15 1.00 2544 2153
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept 47.54 1.71 43.86 50.74 1.00 6544
## CDI.z_age_months 0.07 0.60 -1.09 1.25 1.00 7444
## gender1 -1.99 1.68 -5.33 1.30 1.00 9262
## z_age_months -0.04 0.43 -0.91 0.81 1.00 7060
## method1 0.01 1.32 -2.85 2.90 1.00 6039
## method2 -0.37 1.98 -5.33 3.15 1.00 6133
## IDS_pref -0.66 2.29 -6.80 3.04 1.00 5705
## method1:IDS_pref -0.13 1.99 -4.73 4.07 1.00 7113
## method2:IDS_pref 0.13 2.39 -4.60 5.66 1.00 5917
## CDI.z_age_months:IDS_pref 0.47 1.13 -1.56 3.03 1.00 7947
## Tail_ESS
## Intercept 3678
## CDI.z_age_months 5504
## gender1 5942
## z_age_months 5475
## method1 3608
## method2 2658
## IDS_pref 2326
## method1:IDS_pref 3277
## method2:IDS_pref 2471
## CDI.z_age_months:IDS_pref 4550
##
## Family Specific Parameters:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma 28.65 1.23 26.40 31.20 1.00 7405 6208
##
## Samples were drawn using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).Calculate the marginal log likelihoods for each hypothese using Bridge sampling
marloglik_full_2.nae <- bridge_sampler(lmer.bf_full_2.nae, silent = T)
marloglik_null_2.nae <- bridge_sampler(lmer.bf_null_2.nae, silent = T)Bayes factor
BF_full_2.nae <- bayes_factor(marloglik_full_2.nae, marloglik_null_2.nae)
BF_full_2.nae## Estimated Bayes factor in favor of x1 over x2: 0.810124.7 UK model (main effect of IDS preference)
lmer.bf_full_1.uk <- brm(vocab_nwords ~ CDI.z_age_months + gender +
z_age_months + method + IDS_pref +
(1 | labid),
data = data.uk,
family = gaussian,
prior = prior_uk_1,
iter = 4000,
warmup = 2000,
chains = 4,
save_pars = save_pars(all = T),
control = list(adapt_delta = .95),
seed=213)## Compiling Stan program...## Start sampling##
## SAMPLING FOR MODEL '06bfb59652a0c5924487daf6d3cbc793' NOW (CHAIN 1).
## Chain 1:
## Chain 1: Gradient evaluation took 6.6e-05 seconds
## Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 0.66 seconds.
## Chain 1: Adjust your expectations accordingly!
## Chain 1:
## Chain 1:
## Chain 1: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 1: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 1: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 1: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 1: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 1: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 1: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 1: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 1: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 1: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 1: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 1: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 1:
## Chain 1: Elapsed Time: 0.622256 seconds (Warm-up)
## Chain 1: 0.476119 seconds (Sampling)
## Chain 1: 1.09838 seconds (Total)
## Chain 1:
##
## SAMPLING FOR MODEL '06bfb59652a0c5924487daf6d3cbc793' NOW (CHAIN 2).
## Chain 2:
## Chain 2: Gradient evaluation took 4.1e-05 seconds
## Chain 2: 1000 transitions using 10 leapfrog steps per transition would take 0.41 seconds.
## Chain 2: Adjust your expectations accordingly!
## Chain 2:
## Chain 2:
## Chain 2: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 2: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 2: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 2: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 2: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 2: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 2: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 2: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 2: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 2: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 2: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 2: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 2:
## Chain 2: Elapsed Time: 0.78034 seconds (Warm-up)
## Chain 2: 0.514808 seconds (Sampling)
## Chain 2: 1.29515 seconds (Total)
## Chain 2:
##
## SAMPLING FOR MODEL '06bfb59652a0c5924487daf6d3cbc793' NOW (CHAIN 3).
## Chain 3:
## Chain 3: Gradient evaluation took 7.8e-05 seconds
## Chain 3: 1000 transitions using 10 leapfrog steps per transition would take 0.78 seconds.
## Chain 3: Adjust your expectations accordingly!
## Chain 3:
## Chain 3:
## Chain 3: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 3: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 3: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 3: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 3: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 3: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 3: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 3: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 3: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 3: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 3: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 3: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 3:
## Chain 3: Elapsed Time: 0.679524 seconds (Warm-up)
## Chain 3: 0.907757 seconds (Sampling)
## Chain 3: 1.58728 seconds (Total)
## Chain 3:
##
## SAMPLING FOR MODEL '06bfb59652a0c5924487daf6d3cbc793' NOW (CHAIN 4).
## Chain 4:
## Chain 4: Gradient evaluation took 0.000222 seconds
## Chain 4: 1000 transitions using 10 leapfrog steps per transition would take 2.22 seconds.
## Chain 4: Adjust your expectations accordingly!
## Chain 4:
## Chain 4:
## Chain 4: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 4: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 4: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 4: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 4: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 4: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 4: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 4: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 4: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 4: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 4: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 4: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 4:
## Chain 4: Elapsed Time: 0.740539 seconds (Warm-up)
## Chain 4: 0.838729 seconds (Sampling)
## Chain 4: 1.57927 seconds (Total)
## Chain 4:summary(lmer.bf_full_1.uk)## Family: gaussian
## Links: mu = identity; sigma = identity
## Formula: vocab_nwords ~ CDI.z_age_months + gender + z_age_months + method + IDS_pref + (1 | labid)
## Data: data.uk (Number of observations: 119)
## Samples: 4 chains, each with iter = 4000; warmup = 2000; thin = 1;
## total post-warmup samples = 8000
##
## Group-Level Effects:
## ~labid (Number of levels: 4)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 0.84 0.63 0.04 2.35 1.00 5216 2823
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## Intercept 49.61 3.35 43.36 56.17 1.00 3288 3757
## CDI.z_age_months 22.88 5.82 14.39 34.27 1.00 3814 5617
## gender1 0.83 12.22 -15.11 30.57 1.00 5091 4410
## z_age_months -0.38 1.77 -4.58 2.83 1.00 7242 2973
## method1 -0.41 3.54 -9.07 5.50 1.01 4918 1421
## IDS_pref 0.66 6.25 -7.54 14.21 1.00 2752 968
##
## Family Specific Parameters:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma 155.31 10.18 136.85 176.99 1.00 7865 5540
##
## Samples were drawn using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).lmer.bf_null_1.uk <- brm(vocab_nwords ~ CDI.z_age_months + gender +
z_age_months + method +
(1 | labid),
data = data.uk,
family = gaussian,
prior = prior_uk_null_1,
iter = 4000,
warmup = 2000,
chains = 4,
save_pars = save_pars(all = T),
control = list(adapt_delta = .95),
seed =312)## Compiling Stan program...
## Start sampling##
## SAMPLING FOR MODEL '3d251476a88e0247333a7c190c25e48a' NOW (CHAIN 1).
## Chain 1:
## Chain 1: Gradient evaluation took 6.4e-05 seconds
## Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 0.64 seconds.
## Chain 1: Adjust your expectations accordingly!
## Chain 1:
## Chain 1:
## Chain 1: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 1: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 1: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 1: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 1: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 1: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 1: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 1: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 1: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 1: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 1: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 1: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 1:
## Chain 1: Elapsed Time: 0.494055 seconds (Warm-up)
## Chain 1: 0.475653 seconds (Sampling)
## Chain 1: 0.969708 seconds (Total)
## Chain 1:
##
## SAMPLING FOR MODEL '3d251476a88e0247333a7c190c25e48a' NOW (CHAIN 2).
## Chain 2:
## Chain 2: Gradient evaluation took 5.9e-05 seconds
## Chain 2: 1000 transitions using 10 leapfrog steps per transition would take 0.59 seconds.
## Chain 2: Adjust your expectations accordingly!
## Chain 2:
## Chain 2:
## Chain 2: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 2: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 2: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 2: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 2: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 2: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 2: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 2: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 2: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 2: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 2: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 2: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 2:
## Chain 2: Elapsed Time: 0.581779 seconds (Warm-up)
## Chain 2: 0.51256 seconds (Sampling)
## Chain 2: 1.09434 seconds (Total)
## Chain 2:
##
## SAMPLING FOR MODEL '3d251476a88e0247333a7c190c25e48a' NOW (CHAIN 3).
## Chain 3:
## Chain 3: Gradient evaluation took 7.1e-05 seconds
## Chain 3: 1000 transitions using 10 leapfrog steps per transition would take 0.71 seconds.
## Chain 3: Adjust your expectations accordingly!
## Chain 3:
## Chain 3:
## Chain 3: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 3: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 3: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 3: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 3: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 3: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 3: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 3: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 3: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 3: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 3: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 3: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 3:
## Chain 3: Elapsed Time: 0.520601 seconds (Warm-up)
## Chain 3: 0.405477 seconds (Sampling)
## Chain 3: 0.926078 seconds (Total)
## Chain 3:
##
## SAMPLING FOR MODEL '3d251476a88e0247333a7c190c25e48a' NOW (CHAIN 4).
## Chain 4:
## Chain 4: Gradient evaluation took 8.7e-05 seconds
## Chain 4: 1000 transitions using 10 leapfrog steps per transition would take 0.87 seconds.
## Chain 4: Adjust your expectations accordingly!
## Chain 4:
## Chain 4:
## Chain 4: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 4: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 4: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 4: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 4: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 4: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 4: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 4: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 4: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 4: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 4: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 4: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 4:
## Chain 4: Elapsed Time: 0.553748 seconds (Warm-up)
## Chain 4: 0.495358 seconds (Sampling)
## Chain 4: 1.04911 seconds (Total)
## Chain 4:summary(lmer.bf_null_1.uk)## Family: gaussian
## Links: mu = identity; sigma = identity
## Formula: vocab_nwords ~ CDI.z_age_months + gender + z_age_months + method + (1 | labid)
## Data: data.uk (Number of observations: 119)
## Samples: 4 chains, each with iter = 4000; warmup = 2000; thin = 1;
## total post-warmup samples = 8000
##
## Group-Level Effects:
## ~labid (Number of levels: 4)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 0.84 0.64 0.03 2.37 1.00 5601 3361
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## Intercept 49.67 3.21 43.64 56.03 1.00 4019 3933
## CDI.z_age_months 22.91 5.87 14.36 34.56 1.00 4085 6110
## gender1 0.56 12.09 -15.13 29.92 1.00 4974 4212
## z_age_months -0.38 1.82 -4.81 2.88 1.00 6882 2726
## method1 -0.29 3.16 -7.45 4.94 1.00 5668 2022
##
## Family Specific Parameters:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma 155.74 10.33 137.24 177.67 1.00 7565 4852
##
## Samples were drawn using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).Calculate the marginal log likelihoods for each hypothese using Bridge sampling
marloglik_full_1.uk <- bridge_sampler(lmer.bf_full_1.uk, silent = T)
marloglik_null_1.uk <- bridge_sampler(lmer.bf_null_1.uk, silent = T)Bayes factor
BF_full_1.uk <- bayes_factor(marloglik_full_1.uk, marloglik_null_1.uk)
BF_full_1.uk## Estimated Bayes factor in favor of x1 over x2: 0.940084.8 UK model (interaction between IDS test age and IDS preference)
lmer.bf_full_2.uk <- brm(vocab_nwords ~ CDI.z_age_months + gender +
z_age_months + method + IDS_pref +
IDS_pref:method + IDS_pref:CDI.z_age_months + IDS_pref:z_age_months +
(1 | labid),
data = data.uk,
family = gaussian,
prior = prior_uk_2,
iter = 4000,
warmup = 2000,
chains = 4,
save_pars = save_pars(all = T),
control = list(adapt_delta = .95),
seed=546)## Compiling Stan program...## Start sampling##
## SAMPLING FOR MODEL '5b897001df2f84a1996c21892ffd87d5' NOW (CHAIN 1).
## Chain 1:
## Chain 1: Gradient evaluation took 5.5e-05 seconds
## Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 0.55 seconds.
## Chain 1: Adjust your expectations accordingly!
## Chain 1:
## Chain 1:
## Chain 1: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 1: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 1: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 1: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 1: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 1: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 1: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 1: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 1: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 1: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 1: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 1: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 1:
## Chain 1: Elapsed Time: 0.654481 seconds (Warm-up)
## Chain 1: 0.801309 seconds (Sampling)
## Chain 1: 1.45579 seconds (Total)
## Chain 1:
##
## SAMPLING FOR MODEL '5b897001df2f84a1996c21892ffd87d5' NOW (CHAIN 2).
## Chain 2:
## Chain 2: Gradient evaluation took 7.2e-05 seconds
## Chain 2: 1000 transitions using 10 leapfrog steps per transition would take 0.72 seconds.
## Chain 2: Adjust your expectations accordingly!
## Chain 2:
## Chain 2:
## Chain 2: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 2: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 2: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 2: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 2: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 2: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 2: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 2: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 2: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 2: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 2: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 2: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 2:
## Chain 2: Elapsed Time: 1.04675 seconds (Warm-up)
## Chain 2: 0.558659 seconds (Sampling)
## Chain 2: 1.60541 seconds (Total)
## Chain 2:
##
## SAMPLING FOR MODEL '5b897001df2f84a1996c21892ffd87d5' NOW (CHAIN 3).
## Chain 3:
## Chain 3: Gradient evaluation took 0.000416 seconds
## Chain 3: 1000 transitions using 10 leapfrog steps per transition would take 4.16 seconds.
## Chain 3: Adjust your expectations accordingly!
## Chain 3:
## Chain 3:
## Chain 3: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 3: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 3: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 3: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 3: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 3: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 3: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 3: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 3: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 3: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 3: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 3: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 3:
## Chain 3: Elapsed Time: 1.10452 seconds (Warm-up)
## Chain 3: 0.972033 seconds (Sampling)
## Chain 3: 2.07656 seconds (Total)
## Chain 3:
##
## SAMPLING FOR MODEL '5b897001df2f84a1996c21892ffd87d5' NOW (CHAIN 4).
## Chain 4:
## Chain 4: Gradient evaluation took 6e-05 seconds
## Chain 4: 1000 transitions using 10 leapfrog steps per transition would take 0.6 seconds.
## Chain 4: Adjust your expectations accordingly!
## Chain 4:
## Chain 4:
## Chain 4: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 4: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 4: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 4: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 4: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 4: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 4: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 4: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 4: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 4: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 4: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 4: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 4:
## Chain 4: Elapsed Time: 0.784626 seconds (Warm-up)
## Chain 4: 0.637108 seconds (Sampling)
## Chain 4: 1.42173 seconds (Total)
## Chain 4:summary(lmer.bf_full_2.uk)## Family: gaussian
## Links: mu = identity; sigma = identity
## Formula: vocab_nwords ~ CDI.z_age_months + gender + z_age_months + method + IDS_pref + IDS_pref:method + IDS_pref:CDI.z_age_months + IDS_pref:z_age_months + (1 | labid)
## Data: data.uk (Number of observations: 119)
## Samples: 4 chains, each with iter = 4000; warmup = 2000; thin = 1;
## total post-warmup samples = 8000
##
## Group-Level Effects:
## ~labid (Number of levels: 4)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 0.84 0.64 0.03 2.39 1.00 3969 2427
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept 49.50 3.30 43.15 55.87 1.00 3525
## CDI.z_age_months 22.75 5.82 14.30 34.33 1.00 3826
## gender1 1.05 12.46 -15.70 32.43 1.00 4139
## z_age_months -0.40 1.78 -4.79 2.79 1.00 6214
## method1 -0.25 3.32 -8.11 5.76 1.01 6002
## IDS_pref 0.25 4.67 -8.55 10.43 1.00 5497
## method1:IDS_pref -0.16 4.74 -9.75 8.10 1.00 5314
## CDI.z_age_months:IDS_pref 0.32 3.38 -5.96 8.16 1.00 4842
## z_age_months:IDS_pref 0.10 3.29 -6.72 7.44 1.00 4534
## Tail_ESS
## Intercept 4104
## CDI.z_age_months 6123
## gender1 3799
## z_age_months 3105
## method1 1824
## IDS_pref 1840
## method1:IDS_pref 1657
## CDI.z_age_months:IDS_pref 1953
## z_age_months:IDS_pref 1949
##
## Family Specific Parameters:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma 155.74 10.44 136.75 177.77 1.00 7335 5236
##
## Samples were drawn using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).lmer.bf_null_2.uk <- brm(vocab_nwords ~ CDI.z_age_months + gender +
z_age_months + method + IDS_pref +
IDS_pref:method + IDS_pref:CDI.z_age_months +
(1 | labid),
data = data.uk,
family = gaussian,
prior = prior_uk_null_2,
iter = 4000,
warmup = 2000,
chains = 4,
save_pars = save_pars(all = T),
control = list(adapt_delta = .95),
seed =689)## Compiling Stan program...
## Start sampling##
## SAMPLING FOR MODEL 'fed304567041233af949a5c50d56c33a' NOW (CHAIN 1).
## Chain 1:
## Chain 1: Gradient evaluation took 4.5e-05 seconds
## Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 0.45 seconds.
## Chain 1: Adjust your expectations accordingly!
## Chain 1:
## Chain 1:
## Chain 1: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 1: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 1: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 1: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 1: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 1: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 1: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 1: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 1: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 1: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 1: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 1: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 1:
## Chain 1: Elapsed Time: 0.776728 seconds (Warm-up)
## Chain 1: 0.942291 seconds (Sampling)
## Chain 1: 1.71902 seconds (Total)
## Chain 1:
##
## SAMPLING FOR MODEL 'fed304567041233af949a5c50d56c33a' NOW (CHAIN 2).
## Chain 2:
## Chain 2: Gradient evaluation took 4.4e-05 seconds
## Chain 2: 1000 transitions using 10 leapfrog steps per transition would take 0.44 seconds.
## Chain 2: Adjust your expectations accordingly!
## Chain 2:
## Chain 2:
## Chain 2: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 2: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 2: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 2: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 2: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 2: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 2: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 2: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 2: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 2: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 2: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 2: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 2:
## Chain 2: Elapsed Time: 0.811237 seconds (Warm-up)
## Chain 2: 0.883222 seconds (Sampling)
## Chain 2: 1.69446 seconds (Total)
## Chain 2:
##
## SAMPLING FOR MODEL 'fed304567041233af949a5c50d56c33a' NOW (CHAIN 3).
## Chain 3:
## Chain 3: Gradient evaluation took 5.3e-05 seconds
## Chain 3: 1000 transitions using 10 leapfrog steps per transition would take 0.53 seconds.
## Chain 3: Adjust your expectations accordingly!
## Chain 3:
## Chain 3:
## Chain 3: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 3: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 3: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 3: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 3: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 3: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 3: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 3: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 3: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 3: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 3: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 3: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 3:
## Chain 3: Elapsed Time: 0.744778 seconds (Warm-up)
## Chain 3: 0.957237 seconds (Sampling)
## Chain 3: 1.70202 seconds (Total)
## Chain 3:
##
## SAMPLING FOR MODEL 'fed304567041233af949a5c50d56c33a' NOW (CHAIN 4).
## Chain 4:
## Chain 4: Gradient evaluation took 6.3e-05 seconds
## Chain 4: 1000 transitions using 10 leapfrog steps per transition would take 0.63 seconds.
## Chain 4: Adjust your expectations accordingly!
## Chain 4:
## Chain 4:
## Chain 4: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 4: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 4: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 4: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 4: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 4: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 4: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 4: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 4: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 4: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 4: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 4: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 4:
## Chain 4: Elapsed Time: 1.18464 seconds (Warm-up)
## Chain 4: 1.00483 seconds (Sampling)
## Chain 4: 2.18947 seconds (Total)
## Chain 4:summary(lmer.bf_null_2.uk)## Family: gaussian
## Links: mu = identity; sigma = identity
## Formula: vocab_nwords ~ CDI.z_age_months + gender + z_age_months + method + IDS_pref + IDS_pref:method + IDS_pref:CDI.z_age_months + (1 | labid)
## Data: data.uk (Number of observations: 119)
## Samples: 4 chains, each with iter = 4000; warmup = 2000; thin = 1;
## total post-warmup samples = 8000
##
## Group-Level Effects:
## ~labid (Number of levels: 4)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 0.84 0.64 0.03 2.40 1.00 4572 2694
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept 49.52 3.36 43.12 56.13 1.00 3236
## CDI.z_age_months 22.83 5.80 14.36 34.30 1.00 3651
## gender1 0.84 12.49 -15.51 31.91 1.00 5130
## z_age_months -0.40 1.80 -4.84 2.78 1.00 6954
## method1 -0.28 3.56 -8.65 6.31 1.00 5584
## IDS_pref 0.52 6.32 -7.96 12.25 1.01 4157
## method1:IDS_pref -0.10 5.19 -9.43 8.13 1.01 6863
## CDI.z_age_months:IDS_pref 0.24 2.86 -5.29 6.84 1.00 6084
## Tail_ESS
## Intercept 3398
## CDI.z_age_months 5763
## gender1 4107
## z_age_months 3011
## method1 1793
## IDS_pref 1242
## method1:IDS_pref 2120
## CDI.z_age_months:IDS_pref 2366
##
## Family Specific Parameters:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma 155.84 10.39 136.89 177.51 1.00 8376 5386
##
## Samples were drawn using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).Calculate the marginal log likelihoods for each hypothese using Bridge sampling
marloglik_full_2.uk <- bridge_sampler(lmer.bf_full_2.uk, silent = T)
marloglik_null_2.uk <- bridge_sampler(lmer.bf_null_2.uk, silent = T)Bayes factor
BF_full_2.uk <- bayes_factor(marloglik_full_2.uk, marloglik_null_2.uk)
BF_full_2.uk## Estimated Bayes factor in favor of x1 over x2: 0.917524.9 Combined model (interaction between lanaguage zone and IDS preference)
lmer.bf_full_3.combined <- brm(CDI.prop ~ CDI.z_age_months + language_zone + gender +
z_age_months + method + IDS_pref + IDS_pref:language_zone +
IDS_pref:method + IDS_pref:CDI.z_age_months + IDS_pref:z_age_months +
(1 + CDI.z_age_months | labid),
data = data.uk_nae,
family = gaussian,
prior = prior_combined_full_info,
iter = 4000,
warmup = 2000,
chains = 4,
save_pars = save_pars(all = T),
control = list(adapt_delta = .9999),
seed=100)## Compiling Stan program...## Start sampling##
## SAMPLING FOR MODEL '2ae4a1ae26a91de10302d79941e0a307' NOW (CHAIN 1).
## Chain 1:
## Chain 1: Gradient evaluation took 0.000582 seconds
## Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 5.82 seconds.
## Chain 1: Adjust your expectations accordingly!
## Chain 1:
## Chain 1:
## Chain 1: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 1: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 1: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 1: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 1: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 1: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 1: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 1: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 1: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 1: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 1: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 1: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 1:
## Chain 1: Elapsed Time: 12.3233 seconds (Warm-up)
## Chain 1: 7.96945 seconds (Sampling)
## Chain 1: 20.2927 seconds (Total)
## Chain 1:
##
## SAMPLING FOR MODEL '2ae4a1ae26a91de10302d79941e0a307' NOW (CHAIN 2).
## Chain 2:
## Chain 2: Gradient evaluation took 0.000139 seconds
## Chain 2: 1000 transitions using 10 leapfrog steps per transition would take 1.39 seconds.
## Chain 2: Adjust your expectations accordingly!
## Chain 2:
## Chain 2:
## Chain 2: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 2: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 2: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 2: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 2: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 2: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 2: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 2: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 2: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 2: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 2: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 2: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 2:
## Chain 2: Elapsed Time: 11.6286 seconds (Warm-up)
## Chain 2: 8.87576 seconds (Sampling)
## Chain 2: 20.5044 seconds (Total)
## Chain 2:
##
## SAMPLING FOR MODEL '2ae4a1ae26a91de10302d79941e0a307' NOW (CHAIN 3).
## Chain 3:
## Chain 3: Gradient evaluation took 0.000454 seconds
## Chain 3: 1000 transitions using 10 leapfrog steps per transition would take 4.54 seconds.
## Chain 3: Adjust your expectations accordingly!
## Chain 3:
## Chain 3:
## Chain 3: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 3: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 3: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 3: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 3: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 3: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 3: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 3: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 3: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 3: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 3: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 3: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 3:
## Chain 3: Elapsed Time: 11.724 seconds (Warm-up)
## Chain 3: 8.00645 seconds (Sampling)
## Chain 3: 19.7304 seconds (Total)
## Chain 3:
##
## SAMPLING FOR MODEL '2ae4a1ae26a91de10302d79941e0a307' NOW (CHAIN 4).
## Chain 4:
## Chain 4: Gradient evaluation took 0.000173 seconds
## Chain 4: 1000 transitions using 10 leapfrog steps per transition would take 1.73 seconds.
## Chain 4: Adjust your expectations accordingly!
## Chain 4:
## Chain 4:
## Chain 4: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 4: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 4: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 4: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 4: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 4: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 4: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 4: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 4: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 4: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 4: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 4: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 4:
## Chain 4: Elapsed Time: 12.0478 seconds (Warm-up)
## Chain 4: 13.6773 seconds (Sampling)
## Chain 4: 25.7251 seconds (Total)
## Chain 4:summary(lmer.bf_full_3.combined)## Family: gaussian
## Links: mu = identity; sigma = identity
## Formula: CDI.prop ~ CDI.z_age_months + language_zone + gender + z_age_months + method + IDS_pref + IDS_pref:language_zone + IDS_pref:method + IDS_pref:CDI.z_age_months + IDS_pref:z_age_months + (1 + CDI.z_age_months | labid)
## Data: data.uk_nae (Number of observations: 401)
## Samples: 4 chains, each with iter = 4000; warmup = 2000; thin = 1;
## total post-warmup samples = 8000
##
## Group-Level Effects:
## ~labid (Number of levels: 12)
## Estimate Est.Error l-95% CI u-95% CI Rhat
## sd(Intercept) 0.07 0.04 0.01 0.16 1.00
## sd(CDI.z_age_months) 0.03 0.01 0.02 0.06 1.00
## cor(Intercept,CDI.z_age_months) 0.14 0.38 -0.61 0.78 1.00
## Bulk_ESS Tail_ESS
## sd(Intercept) 1628 1950
## sd(CDI.z_age_months) 3157 5279
## cor(Intercept,CDI.z_age_months) 1855 3392
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept 0.33 0.03 0.28 0.39 1.00 4711
## CDI.z_age_months 0.04 0.01 0.02 0.06 1.00 3982
## language_zone1 0.09 0.04 0.02 0.16 1.00 3801
## gender1 0.04 0.01 0.02 0.06 1.00 9997
## z_age_months -0.00 0.00 -0.01 0.01 1.00 8365
## method1 -0.00 0.03 -0.07 0.07 1.00 4194
## method2 -0.00 0.05 -0.12 0.09 1.00 3421
## IDS_pref -0.01 0.04 -0.08 0.07 1.00 8923
## language_zone1:IDS_pref 0.07 0.05 -0.04 0.17 1.00 7604
## method1:IDS_pref -0.02 0.05 -0.12 0.08 1.00 7668
## method2:IDS_pref -0.05 0.08 -0.21 0.11 1.00 6009
## CDI.z_age_months:IDS_pref 0.01 0.01 -0.01 0.03 1.00 9258
## z_age_months:IDS_pref -0.00 0.01 -0.02 0.02 1.00 8977
## Tail_ESS
## Intercept 3933
## CDI.z_age_months 4494
## language_zone1 3332
## gender1 5754
## z_age_months 6109
## method1 3863
## method2 4142
## IDS_pref 6434
## language_zone1:IDS_pref 5949
## method1:IDS_pref 6349
## method2:IDS_pref 6147
## CDI.z_age_months:IDS_pref 6391
## z_age_months:IDS_pref 6077
##
## Family Specific Parameters:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma 0.21 0.01 0.20 0.23 1.00 8268 5679
##
## Samples were drawn using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).lmer.bf_null_3.combined <- brm(CDI.prop ~ CDI.z_age_months + language_zone + gender +
z_age_months + method + IDS_pref +
IDS_pref:method + IDS_pref:CDI.z_age_months + IDS_pref:z_age_months +
(1 + CDI.z_age_months | labid),
data = data.uk_nae,
family = gaussian,
prior = prior_combined_null_info,
iter = 4000,
warmup = 2000,
chains = 4,
save_pars = save_pars(all = T),
control = list(adapt_delta = .9999),
seed =200)## Compiling Stan program...
## Start sampling##
## SAMPLING FOR MODEL '1e9fb42f69dff2e35aab04bc0c20c833' NOW (CHAIN 1).
## Chain 1:
## Chain 1: Gradient evaluation took 0.000615 seconds
## Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 6.15 seconds.
## Chain 1: Adjust your expectations accordingly!
## Chain 1:
## Chain 1:
## Chain 1: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 1: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 1: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 1: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 1: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 1: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 1: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 1: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 1: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 1: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 1: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 1: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 1:
## Chain 1: Elapsed Time: 13.2027 seconds (Warm-up)
## Chain 1: 8.69344 seconds (Sampling)
## Chain 1: 21.8962 seconds (Total)
## Chain 1:
##
## SAMPLING FOR MODEL '1e9fb42f69dff2e35aab04bc0c20c833' NOW (CHAIN 2).
## Chain 2:
## Chain 2: Gradient evaluation took 0.000373 seconds
## Chain 2: 1000 transitions using 10 leapfrog steps per transition would take 3.73 seconds.
## Chain 2: Adjust your expectations accordingly!
## Chain 2:
## Chain 2:
## Chain 2: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 2: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 2: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 2: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 2: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 2: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 2: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 2: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 2: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 2: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 2: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 2: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 2:
## Chain 2: Elapsed Time: 12.3354 seconds (Warm-up)
## Chain 2: 13.267 seconds (Sampling)
## Chain 2: 25.6023 seconds (Total)
## Chain 2:
##
## SAMPLING FOR MODEL '1e9fb42f69dff2e35aab04bc0c20c833' NOW (CHAIN 3).
## Chain 3:
## Chain 3: Gradient evaluation took 0.000189 seconds
## Chain 3: 1000 transitions using 10 leapfrog steps per transition would take 1.89 seconds.
## Chain 3: Adjust your expectations accordingly!
## Chain 3:
## Chain 3:
## Chain 3: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 3: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 3: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 3: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 3: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 3: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 3: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 3: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 3: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 3: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 3: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 3: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 3:
## Chain 3: Elapsed Time: 12.4057 seconds (Warm-up)
## Chain 3: 7.8485 seconds (Sampling)
## Chain 3: 20.2542 seconds (Total)
## Chain 3:
##
## SAMPLING FOR MODEL '1e9fb42f69dff2e35aab04bc0c20c833' NOW (CHAIN 4).
## Chain 4:
## Chain 4: Gradient evaluation took 0.000322 seconds
## Chain 4: 1000 transitions using 10 leapfrog steps per transition would take 3.22 seconds.
## Chain 4: Adjust your expectations accordingly!
## Chain 4:
## Chain 4:
## Chain 4: Iteration: 1 / 4000 [ 0%] (Warmup)
## Chain 4: Iteration: 400 / 4000 [ 10%] (Warmup)
## Chain 4: Iteration: 800 / 4000 [ 20%] (Warmup)
## Chain 4: Iteration: 1200 / 4000 [ 30%] (Warmup)
## Chain 4: Iteration: 1600 / 4000 [ 40%] (Warmup)
## Chain 4: Iteration: 2000 / 4000 [ 50%] (Warmup)
## Chain 4: Iteration: 2001 / 4000 [ 50%] (Sampling)
## Chain 4: Iteration: 2400 / 4000 [ 60%] (Sampling)
## Chain 4: Iteration: 2800 / 4000 [ 70%] (Sampling)
## Chain 4: Iteration: 3200 / 4000 [ 80%] (Sampling)
## Chain 4: Iteration: 3600 / 4000 [ 90%] (Sampling)
## Chain 4: Iteration: 4000 / 4000 [100%] (Sampling)
## Chain 4:
## Chain 4: Elapsed Time: 11.9008 seconds (Warm-up)
## Chain 4: 8.74046 seconds (Sampling)
## Chain 4: 20.6412 seconds (Total)
## Chain 4:summary(lmer.bf_null_3.combined)## Family: gaussian
## Links: mu = identity; sigma = identity
## Formula: CDI.prop ~ CDI.z_age_months + language_zone + gender + z_age_months + method + IDS_pref + IDS_pref:method + IDS_pref:CDI.z_age_months + IDS_pref:z_age_months + (1 + CDI.z_age_months | labid)
## Data: data.uk_nae (Number of observations: 401)
## Samples: 4 chains, each with iter = 4000; warmup = 2000; thin = 1;
## total post-warmup samples = 8000
##
## Group-Level Effects:
## ~labid (Number of levels: 12)
## Estimate Est.Error l-95% CI u-95% CI Rhat
## sd(Intercept) 0.06 0.04 0.01 0.15 1.00
## sd(CDI.z_age_months) 0.03 0.01 0.02 0.06 1.00
## cor(Intercept,CDI.z_age_months) 0.11 0.39 -0.67 0.77 1.00
## Bulk_ESS Tail_ESS
## sd(Intercept) 1755 3160
## sd(CDI.z_age_months) 3490 4973
## cor(Intercept,CDI.z_age_months) 1652 3129
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept 0.33 0.03 0.28 0.39 1.00 4857
## CDI.z_age_months 0.04 0.01 0.02 0.06 1.00 3471
## language_zone1 0.09 0.03 0.02 0.16 1.00 3922
## gender1 0.04 0.01 0.02 0.06 1.00 10010
## z_age_months -0.00 0.00 -0.01 0.01 1.00 8804
## method1 -0.00 0.03 -0.07 0.06 1.00 4326
## method2 -0.00 0.05 -0.11 0.09 1.00 3617
## IDS_pref -0.01 0.04 -0.08 0.07 1.00 7461
## method1:IDS_pref -0.02 0.05 -0.12 0.09 1.00 6112
## method2:IDS_pref 0.01 0.06 -0.11 0.14 1.00 6422
## CDI.z_age_months:IDS_pref 0.01 0.01 -0.01 0.03 1.00 9622
## z_age_months:IDS_pref -0.01 0.01 -0.03 0.02 1.00 7741
## Tail_ESS
## Intercept 3521
## CDI.z_age_months 4281
## language_zone1 3649
## gender1 6022
## z_age_months 6358
## method1 4570
## method2 4296
## IDS_pref 6270
## method1:IDS_pref 6007
## method2:IDS_pref 5734
## CDI.z_age_months:IDS_pref 5430
## z_age_months:IDS_pref 5676
##
## Family Specific Parameters:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma 0.21 0.01 0.20 0.23 1.00 8717 5608
##
## Samples were drawn using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).Calculate the marginal log likelihoods for each hypothese using Bridge sampling
marloglik_full_3.combined <- bridge_sampler(lmer.bf_full_3.combined, silent = T)
marloglik_null_3.combined <- bridge_sampler(lmer.bf_null_3.combined, silent = T)Bayes factor
BF_full_3.combined <- bayes_factor(marloglik_full_3.combined, marloglik_null_3.combined)
BF_full_3.combined## Estimated Bayes factor in favor of x1 over x2: 0.09910```